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ABSTRACT 


Frequently, in struetural system identifieation (model updating or damage detection), the 
available set of data is incomplete, both spatially and in modal content. This 
incompleteness leads to the solution of underdetermined linear systems. In order to 
improve the identification, additional independent measured data must be found. In the 
past, it has been shown that such data can be easily obtained from the application of 
Artificial Boundary Conditions (ABC), imposed on both the baseline FE models and the 
measured frequency response data. This can be accomplished without any physical 
modifications to the experiment and, hence, no additional expense on different systems, 
or more than once, in order to get the modal data needed for the analysis. In this thesis, 
the procedure of sensitivity-based structural system identification, using ABCs, and 
enhanced by parameter grouping/clustering, is examined. It is shown that the optimal 
sensitivity matrix is a square and diagonal dominant one, which can be used with quite 
accurate results both for localization of parameter errors, and the determination of the 
magnitude of parameter error. The numerous ABC configurations available, even from a 
small measured data set, allow an optimal sensitivity matrix to be found for many 
parameters. These concepts are demonstrated using simulated measurements along with 
finite element models. 
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I. 


INTRODUCTION 


A. MODEL UPDATING AND DAMAGE DETECTION 

Finite element analysis (FE) is a computational tool that is widely used in cases 
where a closed-form solution to a partial differential equation is not available. In 
structural dynamics, FE models are used to predict the dynamic response of a system to 
external stimuli, such as dynamic excitations, or to assess the effect of parameter 
changes. The term “parameter” stands for all the individual physical and geometric data 
that define the structure, such as geometry, elastic moduli and element masses/densities. 

The use of EE models aims at achieving higher accuracy in terms of response 
prediction; however, doing so may not always be possible for several reasons. Some of 
them are usually connected to uncertainties regarding the governing equation of the 
system, to assumptions made during the derivation, to inaccurate physical boundary 
conditions, to inexact characterization of joints and connections, etc. All these 
discrepancies are revealed and attempted to be overcome, when a necessary procedure for 
the validation of the EE model takes place. This procedure is called EEM updating. It 
involves measuring structure response, identifying modal data, comparing the modal 
parameters to those extracted from the EE model and finally, modifying appropriately the 
EE model to match the system’s design requirements. EEM updating is of vital 
importance as it verifies that the model can be used with acceptable reliability for a 
number of research fields in modal analysis, such as damage detection and error 
localization. 

B. BACKGROUND 

Over the last decades, a large volume of literature on structural identification, 
including both model updating and error localization has been produced. Eriswell and 
Mottershead in their survey [1], extensively discuss methods that deal with various 
aspects of this field. Most of them are eigenvalue sensitivity-based and are thoroughly 
presented and analyzed in [2]. Each developed method has, as a main objective, to 
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overcome problems that arise during the procedure, which in general, converges to the 
solving of a linear system. The most common issue with the solution is the ill- 
conditioning of these systems. This ill-condition results from the lack of input data, 
which frequently leads to a severely underdetermined system of equations. The lack of 
input data is mainly caused either by the geometric limitation of measurement locations 
on a structure (spatial incompleteness) or by the bounds of the available bandwidth 
during modal testing (modal incompleteness). It is often observed that ill-conditioning 
may be a concern due to dependencies of the sensitivity matrix, or even due to the 
physical geometry of the structure, in cases where symmetry is present. 

Dealing with this problem, Lallement and Piranda in [3] are the first to introduce 
the development of a rational procedure for parameter selection—^which appears to be the 
main core in many updating methods—to facilitate the improvement of the conditioning 
of the equations. This has to do with an iterative procedure to find an appropriate subset 
of parameters that produce an optimal solution commonly known as forward selection. 
Then, the importance of vector subspaces and their wide application range in modal 
analysis is revealed [4], and in [5] the concept of selecting parameters groups using 
angles between subspaces is implemented, also incorporating side constraints to the 
equations to be solved. The methods of subset selection are outlined and extended to the 
selection of groups of parameters. Reference [6] is a representative example of the 
selection and updating of parameters in an actual application. Another approach to the 
same objective, parameter selection, is offered in [7] where an automated procedure to 
guide the parameter selection is suggested based on simple observations. The idea of 
grouping element sensitivities together, under the condition of selecting as few 
parameters as possible that comply with the above stated observations, is presented. 
These observations deal with the fact that the number of updating parameters can be 
reduced, but only by sacrificing the total sensitivity, which is a concept that will be 
covered and used in Chapter IV of this thesis. References [8], [9] and [10] address this 
concept in a different manner, using the principles of cluster analysis, along with optimal 
subset selection, and show how these methods can be successfully used in a number of 
applications. 
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Another path towards the improving the ill-eonditioning aims at the enrichment of 
the available data, in order to decrease incompleteness, or be less underdetermined. Of 
course, care has to be exercised in preventing the addition of redundant data and in 
accordance with feasibility conditions. A first step of doing so without requiring any 
modifications to the physical structure is made by Gordis [11], where the concept of 
artificial constraints is initially introduced. The term artificial refers to the fact that a 
constraint is applied computationally to the measured data, as well as to the FE model. 
Another similar reference is made by Rade, Lallement and Silva in [12], stated here as 
Modified Configurations, describing related ways different artificial systems are used to 
enlarge the experimental data. In this work, the criterion established is the minimization 
of the condition number of the sensitivity matrix, or in other words, the maximization of 
the degree of linear independence between the columns and the rows of the sensitivity 
matrix. A theoretical investigation of this concept is made in [13], where the term ABC, 
standing for Artificial Boundary Condition, is first introduced. It proves that, with no 
additional testing, all the natural frequencies corresponding to ABC systems and the 
needed additional data can be obtained from any square FRF matrix measured in a modal 
test. In addition, it shows that structural sensitivities can be generated from both the 
baseline FE model and from the ABC systems, and combined in a global sensitivity 
matrix which can reduce or eliminate conditioning problems. A genetic algorithm (GA)- 
based methodology, to make effective use of the ABC frequencies for FE model 
updating, is presented in [14] by Fu and Tu, where a binary code, GA, is proposed for the 
selection of the desired artificial boundary conditions. The recent literature related to 
ABCs in general is completed with [15], where the practical aspects of the 
implementation of ABCs in an actual measurement environment is investigated and the 
feasibility of using only the lowest few ABC frequencies in conjunction with a variety of 
ABC configurations is examined. 

C. REFERENCE TO PRESENT WORK 

In this thesis, a sensitivity-based FE analysis is performed, with the goal of 
eliminating the discrepancies mentioned above. First of all, the properties of several 


3 



Matlab algorithms are examined in order to choose the one that best matches our 
requirements for linear system solution. These algorithms are invoked with the 
commands “pinv” and ‘\” and, eventually, involve SVD and QR decomposition with 
column pivoting, respectively. In addition, the correlation between the sensitivity matrix 
rows and modal strain energy is depicted. Then, the degree of the data incompleteness, 
which is related to the number of available modes, is connected to the qualitative result of 
the linear system solution. This leads further to the conclusion that there is a necessary 
but not sufficient condition for detecting potential damage everywhere in a structure, 
using only one sensitivity matrix for all the cases. The condition is that one must provide 
as many modes as the adjustable design parameters of the structure, i.e., building a square 
sensitivity matrix. Toward this end, several methods of data enrichment are taken into 
account, like eigenvector sensitivities and clustering/grouping parameters, with emphasis 
on the use of ABCs, since they are considered the most effective tool for all the 
aforementioned reasons. Keeping in mind that even with a potentially large number of 
ABC configurations available, a sufficient amount of additional data from the ABCs may 
not be sufficient, such as in very large and complex structures; in these cases, ABC 
configurations can be combined with grouped parameters into a global sensitivity matrix. 
Finally, it is shown that the optimal configurations result in a diagonal dominant global 
sensitivity matrix, which once built, will meet necessary and sufficient conditions for 
damage detection. 
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II. THEORY 


Beginning from the governing equation that describes a dynamic system under an 
excitation, we have: 

[M]{x} + [C]{x} + [ii:]{x} = {F(t)} (2.1) 

where this system consists of mass, damping and stiffness parameters, expressed as the 
symmetric mass, damping and stiffness matrices shown above, of size n x n, where n is 
the number of degrees of freedom (DOF) in the model. The vector represents the 

excitation force of size n xl. 

If the damping [c] is neglected. Equation (2.1) becomes: 

[M]{x} + [^]{x} = {E(0} (2.2) 


Expressing in a more explicit way each matrix shown in Equation (2.2), results in 
the following: 



(2.3) 


Equation (2.3) describes a system of n degrees of freedom, which approximates a 
continuous system with an infinite number of degrees of freedom. An actual modal test 
involves the measurement of only a small percentage of a structure’s DOEs. During this 
procedure, these measured coordinates are considered the analytical set or (ASET), and 
they are the locations where the measuring instruments (accelerometers) are positioned. 
The remaining unmeasured DOE are considered the omitted coordinate set or (OSET), 
and in theory, contain an infinite number of degrees of freedom. 
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A. 


ANALYTICAL AND OMITTED COORDINATE SETS 


Going from the time domain to the frequeney domain and getting the steady state 
response of the above described system, we get for the whole set of n DOF: 
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(2.4) 


where Q is the frequency of the harmonic excitation. Making use of the analytical 
(measured) and the omitted (unmeasured) set of coordinates, that is the ASET and OSET 
coordinates. Equation (2.4) takes the following form: 


K 

K 


aa 

oa 





(2.5) 


Without loss of any generality, we can assume that there is no excitation on the omitted 
coordinates, that is {/^} = {0} , and get: 
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( 2 . 6 ) 


The above matrix form represents a system of two equations. Erom those equations, there 
is an obvious solution for the OSET coordinates in terms of ASET coordinates through: 

]{x„ } + [if „„ ]{x„} - O^ [[m „„ ]{x„} + [m „„ ]{x„}] = 0 (2.7) 

Rearranging and solving for|x^} , Equation (2.7) becomes: 

K} - [/ ]■' [- ]k} (2.8) 

By definition, the term [^7 - J is equal to: 
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(2.9) 


\1-Ql^K^M T' = 




:Adj\j 




where the DetTz-Q^^fir^^ indieates the determinant and Ad)T/-Q^Zr^„ 

indieates the adjoint matrix. From Equation (2.9), it is obvious that there are some 
frequency values for which the inverse does not exist, and these are the ones that satisfy: 

( 2 . 10 ) 


It is easily noticeable that the frequencies that satisfy the above equation are the 
eigenvalues of the system defined by the [ZT^^jand \M^J\ matrices, hence, they 
correspond to the natural frequencies of the OSET system. 

Realizing that the eigensolution of a dynamic system is defined by the 
unrestrained DOE of the corresponding system, the OSET coordinates are unrestrained 
and the ASET coordinates are considered constrained to the ground. 

B. REDUCED ORDER MODEL 


The fact that a finite number of transducers are used for a modal test, in contrast 
to the infinite number of DOE of the real system, defines a reduced order model. This 
model’s impedance is nonlinearly dependent on the impedance of the full order model 
[13]. Considering this exact model, in terms of frequency response function (ERE) we 
have: 


[H] 
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( 2 . 11 ) 


The experimental set of data measured in a real test is carried in the partition 
of the above full ERE. In other words, [ZZ^^] represents a structural dynamic model, 

reduced via exact dynamic reduction [11]. Using the algebraic relation of ERE to its 
impedance, we have: 
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( 2 . 12 ) 
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which, after decomposing the set of equations and solving for ], leads to: 




-Z,, 



(2.13) 


Writing Equation (2.6) in terms of impedance, we get: 



(2.14) 


Rearranging (2.14), the exact dynamic relationship between the OSET coordinates and 
the ASET coordinates gives: 





(2.15) 


Combining Equation (2.14) and Equation (2.15), and solving the first equation of the 
matrix formula, the dynamic reduced model results: 

{/| = rz -Z Z-‘Z I (2.16) 

Paralleling Equations (2.13) and (2.16), it is found that the ERE extracted from an 
actual modal test, i.e., the measured partition of an infinite-order ERE matrix, is 
equivalent to the inverse of a dynamically reduced impedance matrix. In addition, the 
term ' seems to be critical for the definition of[//^^], and since every element in the 

' is singular at the natural frequencies of the OSET, Equation (2.13) shows that the 
elements of ' will be singular at the OSET natural frequencies [13]. 
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From what is proved above, it is obvious that the FE model can be manipulated in 
the same manner. This means that the FE model is used to calculate a full n x n ERE, and 
appropriate boundary conditions can be applied corresponding to the coordinates 
measured in a test. 

C. FREQUENCY RESPONSE FUNCTION 

The Frequency Response Function FRF is actually a ratio defined by the response 
of a system due to a harmonic excitation. It contains all the amplitude and phase 
information for the response of a specific DOF excited by an input force at a specific, and 
perhaps different, DOF. The FRF is constructed by measuring both the response and the 
excitation simultaneously. This data is transformed from the time domain to the 
frequency domain using FFT algorithm. Doing this for all the DOF that are set to be 
measured for certain external excitations, the FRF matrix [//(Q)] can be built. As 

stated before, the inverse of this matrix is known as the impedance matrix [Z (O)] : 

[Z(Q)] = [H(Q)] ‘ (2.17) 


where: 


[Z(Q)] = [A-Q"M + ;QC] 


(2.18) 


In terms of an n DOF system and introducing impedance. Equation (2.4) can be written 
as: 



_1 

\ 


7i' 

. . . 

N 

•^nn _ 

\ 

> = < 

... 


(2.19) 


Fetting{jc} = [®]{<?} , where O is the mass normalized mode shapes of the system and 
<7 represents the generalized coordinate set, substituting in (2.19): 

[zl<i>fej = {/} 

9 


( 2 . 20 ) 



Then premultiplying by O]^ and combining Equations (2.18) and (2.20); 


[d>f [ii:][cI)]-Q^ [cbf [M][d>]+ 7 Q[cl)f [C][cl)] = [cbf {f} (2.21) 


Since O is mass normalized[O]^ [M][0] = 1, and assuming proportional damping 

[C]^a^[K] + P*[M\. 


-Q^ + 274 nu>,]{^} = [of {F} (2.22) 


In Equation 2.22, co. is the i*’" natural frequency, Q. is the forcing frequency and 
is the systems damping ratio. The term +2ji^cOi'^is known as the modal 

impedance matrix and is diagonal. Inverting the modal impedance matrix and moving 
back to physical coordinates by premultiplying by [O], the frequency response is: 


= [o] 


G)f -fP + 2 


[0f{F} 


(2.23) 


Erom Equation (2.23), the ERE matrix [// (Q)] can be defined as: 


[//(Q)]=[o] 


cof +2 


[Of 


(2.24) 


The ERE matrix can be rewritten in summation form: 

{o^fo* 


modei' 


[//(«)]= £ 


2^ CO), ~ + 2 y ^52 cOf, 


(2.25) 


and descriptively for any element in terms of H.j : 
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(2.26) 


[».(»>]-1 . tSL 

k=\ cOj^— Cl +2j(^ClcOi^ 

It is obvious from Equation (2.26) that the FRF contains all the modal data for the 
system. 


1. Driving Point and Transfer Point Functions 


The exact FRF matrix, as mentioned, contains all the modal data of the system 
and thus the resonance and anti-resonance frequencies, which are of particular interest in 
modal analysis. It should be noted here that, depending on the relationship, in terms of 
locality, of the i'*' DOF where the response is measured with the f' DOF where the force 
is applied, the FRF is referred to as the Driving Point Function for the case that these two 
DOF are identical and Transfer Function for all the other cases. 

What follows is an example of FRF derivation, covering both Driving Point and 
Transfer Function, for a simple two DOF system with no damping [16]. 



Figure 1. 2-DOF Example 


The properties of this system are summarized in the stiffness and mass matrices 


above: 




K,+K, -K, 

-K, K,+K, 


and [M ] = 


Mj 0 

0 M, 


Fet, for this example, Mj =M 2 =1.0 and =0.4 and =0.8. This yields the 

following results: 
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Mode shapes: {d>'| = 


-0.7071| 

-0.70711 


and = 


[-0.7071] 

0.7071 


, , f0.6325] 

Natural frequeneies: [coirad / see)| = < ^ f 


mode^ 


Using these values in Equation (2.26), Hy{D.)^ = ^ 




(U,-Q +27^0®, 


- results in: 




{of}{o(}'^ _ {-0.7071}{-0.7071}'' {-0.7071}{-0.7071}' 


co^-D.^ coi-Q^ 


0.4-Q" 


1.2-Q" 


(2.27) 


for the Driving Point funetion and for the Transfer funetion in: 


Th ctI {-0.7071}{-0.707lf j-0.707lH0.7071}' 

^ 0 ) 1 - 0 ^ 0.4-Q" 1.2-Q" 


(2.28) 


What is evitable from Equation (2.27) for the driving point is that for frequeneies 
below and above the natural frequencies of the system, the two terms will have the same 
sign, while between them, they will have the opposite sign. This means that in the first 
case, the two terms are additive, where in the second they are subtractive. Eurthermore, 
this effect creates resonances at forcing frequencies equal to the natural frequencies and 
for the particular values between them, the two terms of the equation may result in zero 
value, thus points of anti-resonance are created. Therefore, anti-resonance follows every 
resonance, without exception, as seen in the top portion of Eigure 2. 

On the other hand, there is no reason to expect something similar for the Transfer 
Eunction, because obviously the two terms of Equation (2.28) will sometimes be additive 
and sometimes be subtractive, as the product of the two eigenvectors—called modal 
constant—will be sometimes positive and sometimes negative. Thus, we expect Transfer 
ERF measurements to show a mixture of anti-resonances and minima. It can be shown 
that, in general, the further apart the two points in question are, the more likely the two 
eigenvector elements are to alternate in sign as one progresses through the modes [17]. 
The above conclusions are shown in the lower portion of Figure 2. 
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Driving Point FRF, H11 



Transfer FRF, H12 



Figure 2. 2-DOF Frequency Response Function. 

The principles demonstrated in the above example can be applied to any number 
of DOF. This concept of the existence of anti-resonance between any two modes of a 
driving point FRF is of great significance in the use of Artificial Boundary Conditions. 

D. ARTIFICIAL BOUNDARY CONDITIONS 

As discussed previously, the generation of additional data for a spatially 

incomplete FE model analysis can be pursued via Artificial Boundary Conditions, or 

ABC. This way, a number of ABC can be applied to an FE model, thus constructing 

several models differing only in boundary conditions. Hence, a new model configuration 

is provided, yet only one set of measured data is required. Applying the driving point 

function concept shows that the ABC frequencies of a certain ASET coordinate are 

defined by the corresponding driving point anti-resonances. In other words, these ABCs 

are the boundary conditions that define the OSET, from which it is shown that additional 

frequency information about the model is provided. The following examples describe 
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how anti-resonance frequencies are identical to the natural frequencies corresponding to 
the system with the associate driving point constrained to the ground [13]. 


1. ABC Frequencies for a 2-DOF System 

Consider the system shown in Figure 3, excluding damping: 



Figure 3. 2-DOF System 


From Equation (2.26), the driving point FRF is 




modes 


k=l 


0)1 


(2.29) 


where d>, is the mass normalized mode shape element, co/t is the A:th natural frequency, 
and Q. is the forcing frequency. 

Solving for the frequency of the anti-resonance of //jj(Q) yields: 


n: 


RI^co2 + RfiCOi 

Rl+Rn 


(2.30) 


where the modal residue is given by 7?* = cp^cp^. 

Let Ml = M2 = 1.0 and K1 = K2 = 1.0. Then the stiffness and mass matrices of 
the system are: 



-Kj 


■ 1 -\ 


[M, 

0 " 


"1 0“ 

1 I 

1 

— 


and \M ] = 

1 


— 


_ -K2 

K, + K2_ 


-1 2 

0 



0 1 


Solving the eigenproblem related to them yields to: 
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Natural frequencies: [co{rad / sec)} = < 
and the single anti-resonance 


0.6181 

1 . 618 J 

= V 2 rad/sec, as shown in Figure 4. 


O 

o 

o 


Driving Point FRF, Fi11 



rad/sec 


Figure 4. Driving Point FRF, HI 1 
Now, consider the same system pinned at DOF 1, as shown in Figure 5: 
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Figure 5. DOF #1 Constrained to Ground 

This system’s natural frequency is co = y[2 rad/sec and identically equal to the 
anti-resonance of the unconstrained system of Figure 3. The same result for the 
constrained system can be found by making use of Equation (2.13), that is, by calculating 
the inverse of the ASET ERE partition that is constrained//(Q), which in this case is 

the inverse of//jj(Q). The correlation between the driving point anti-resonance and the 
natural frequency of the system with the driving point constrained to ground is shown in 
Eigure 6 [13]. 


Plot A. Driving Point FRF, H11 



rad/sec 

Plot B. Z=inv{H11) of the ABC System 



.50 1 _I_[_^_I_I_[_1_ I _I_I 

0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 

rad/sec 


Figure 6. Plot A. Driving Point HI 1(D) of System 1, Plot B. [Haa(D)]-l of System 2 
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2. 


ABC Frequencies for a Free-Free Beam 


The model depicted in Figure 7 consists of ten beam elements and each element 
contains two nodes with two DOF each, the odd-numbered DOF being translational, and 
the even numbered DOF being rotational. Response transducers are set up at the 
translational DOF; 1, 5, 9, 13, 17, and 21 as follows: 



Figure 7. Transducers Located at DOF 1, 5, 9, 13, 17, and 21 

The location of the response accelerators represent the ASET, defined as [1, 5, 9, 
13, 17, 21], If it is assumed that the excitation is applied at each of the ASET DOE, then a 
square (6 x 6) ERE matrix will result from the testing. Then the impedance matrix 
is calculated in the bandwidth 0-800 Hz. Eigure 8 shows the driving point ERF 
of the system. 



Frequency (Hz) 


Figure 8. Driving Point FRF ASET DOE: 1 5 9 13 17 21. (Erom [13]) 
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Now, consider the same system with all the ASET DOF constrained to the 
ground, as shown in Figure 9. 



i I 1 I i I 

Figure 9. ABC Configuration ASET [1 5 9 13 17 21] (Measured Coordinates Restrained 

to Ground) 

By calculating the impedance matrix and plotting the of the matrix in the 

same bandwidth of 0-800 Hz, Figure 10 is generated; 



Figure 10. Impedance FRF for ASET DOF: 1 5 9 13 17 21 (From [13]) 

By deriving the ABC frequencies of the system in Figure 9, it is obvious that they 
correspond exactly to the resonances of the impedance FRF for ASET coordinates the 
DOF, which are constrained to the ground. 
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E. SENSITIVITY BASED FEM UPDATING 


Sensitivity analysis in structural dynamics is a tool used for predicting the effect 
of small changes in design parameters on dynamic response characteristics. It involves 
finding how a system alters its dynamic response when some parameters that define its 
properties are subject to a change. From a mathematical point of view, it implies 
differentiation, which in the present analysis, is carried out through finite difference 
calculations. The parameters of the system to be changed, stated as design variables, are 
in this case, only stiffness parameters, whereas other parameters, like mass or density can 
be taken into consideration. The output of the system to be examined using sensitivities, 
due to an input parameter change, can be one or both of the following: natural 
frequencies, eigenvalue sensitivity and mode shapes, eigenvector sensitivity. Another 
aspect of modal data analysis via sensitivity is the Frequency Response Function 
Sensitivity Analysis, which is not covered in this thesis. 

Recalling the ABC concept, it is shown that ABC configurations can be used 
along with the baseline system into a combined sensitivity matrix, for model updating 
[13]. The key part of this is that ABC frequencies can easily be found, as stated 
previously, and since the design variables are common for both the ABC and the baseline 
system, they can be used to reliably enlarge the modal data for the analysis. 

The equation that describes what is mentioned above about sensitivity based 
updating is: 

{A®"}=[r]{ADV} (2.31) 


where is a vector containing the difference between the experimental and the 

analytical (FE baseline) natural frequencies ~^a} ■ ^he {ADV} term is the vector 
containing the design variables, and [T] is the sensitivity matrix. 

Each column of the sensitivity matrix represents an element of the model, 
whereas each row represents a mode of the model. 
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Each ABC system defines additional rows to Equation (2.31). It is obvious from 
this point that even with low modes, in terms of frequeney speetrum, that is with the first 
few modes of baseline and ABC systems, a more determined sensitivity matrix ean be 
built. Representing the baseline system quantities with the superseript “0” and the ABC 
systems with the superscripts from 1 to “k,” where k represents the number of ABC used 
in the system, the following equation shows how this is eompromised [13]; 


A® 

Ao 


2 ( 0 ) 

BASE 

2 ( 1 ) 

'ABC\ 


^7’(0) A 

^ BASE 
T’)!) 

^ ABC\ 



(2.32) 


Ao 


2{k) 

ABCk 


rp(k) 

ABCk J 


dV. 


1. Eigenvalue Sensitivity Matrix Derivation. 

The starting point of the derivation is the eigenvalue problem of an n DOE 

system; 




(2.33) 


Differentiating with respeet to the design variable {D V}; 


ADV 


, AM 

— /it- 

ADV 


Ah 

ADV 


{^i} + [K - AiM] 


AO, 

ADV] 


(O) 


(2.34) 


Expanding Equation (2.34) and premultiplying by leads to; 




ADV 


{<!>,}-/li {<5,}^ 


AM 

ADV 


O, 


ADV 


[M]{d),} +{(5,}^ [^_ = {0} (2.35) 


By invoking the matrix identity property [^]{^} = [^]{^} ’ 

left hand side ean be replaeed by; 


AO) Y 
ADV J 


[7f-/l,M]{0,} = 0 


(2.36) 
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Since = |0} and recalling the mass normalization of the mode shapes, 

where overall equation is redueed to: 


O, 


AK 

ADV 




AM 

ADV 


I®'}- 


AAi 

ADV 


= { 0 } 


(2.37) 


Bringing the mass and stiffness portion to the right hand side 


" AAi ' 


■ AK 

, AM 




Ai 

_ADV _ 

ADV 

ADV 


{*,) 


(2.38) 


As stated previously, the present analysis only takes into account changes in the 

stiffness parameters. Therefore, the seeond term Ai in the square brackets of the 

ADV 

right hand side ean be negleeted as AM = 0, yielding to: 

{Oi} where [A^] = [&-^a] (2.39) 


AAi 

ADV 


= 10 , 


AK 

ADV 


and finally: 




AK 

ADV 


{*,) 


(2.40) 


It has been shown [18] that when ABC are applied, the model tends to concentrate 
high stiffness sensitivity values around the partieular DOF where the ABC is loeated. 
This concept, along with the fact that the baseline model is not likely to prediet ehanges 
in natural frequencies in regions with low sensitivity values, is of great significance in the 
analysis of the following chapters. 

2. Eigenvector Sensitivity Matrix Derivation 

The ealculation of the eigenveetor sensitivity is a less straightforward problem 
than that of the eigenvalue sensitivity. The foundation for the ealculation was laid 
primarily by Fox and Kapoor [19] and a simplified method, whieh is the basis of the 
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NASTRAN computation scheme, is provided later by Nelson [20]. A survey of several 
eontributions to the literature regarding the eigenveetor sensitivity ealeulation is available 
in [21], 

The derivation that follows is based on Blelloeh’s report [22], After having 
ealeulated the eigenvalue sensitivity from Equation (2.39), returning to Equation (2.34) 
and solving for the eigenveetor derivative, again for no mass variation yields to: 


[k-xm' 


^(j)i ' 

M)V J 


A/f 

^DV 



ADV 


Ml 


(2.41) 


Assume at this point that all the eigenveetors are available for the system under 
eonsideration. In this ease, the eigenveetor derivative ean be expressed as a linear 
eombination of all the eigenveetors: 

{^} = [®]W P-42) 


where [d)] is the matrix eontaining all the eigenveetors and |c} is the veetor of the 
eoeffieients to be determined. Substituting to Equation (2.41) and premultiplying by 
[or generates: 


[of [7f-yl;M][0]{c} 



AK 

ADV 


AAi 

ADV 


M 


Ml 


(2.43) 


The expression that premultiplies |c} is of the form: 

[O]^ [K - TiM ] [O] = [A - 2 , 1 ] 


(2.44) 


where A is the diagonal matrix of eigenvalues. The above implies that, for distinet 
eigenvalues, every element of the eoeffieient veetor |c} exeept for the ith element (cj is 
uniquely determined as follows: 
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(2.45) 


c 


k 



AK AAi 

- M 

ADV ADV 

K ~ \ 



At this point, it should be noted that there is no reason to believe that the 
eigenveetor derivative will prediet the ehange in the eigenveetor in models with elosely 
spaeed eigenvalues, sinee from what can be seen from Equation (2.45) numerical ill- 
conditioning will occur in these cases. 

The eigenvector derivative can now be written as: 

■ = E {A) + = {r} + ) (2.46) 


In this expression, the vector } is determined uniquely. The undetermined coefficient 
c, will be found on the basis of differentiating the mass normalization formula for the 
eigenvectors: 


to have: 






dM 

dDV 


M} = o 


(2.47) 


Again, for no mass variation and substituting Equation (2.46) into Equation (2.47): 

(2.48) 

The key to the whole method that Nelson [20] indicated is that there is no need to 
calculate the unique|l/.}. In fact, any vector that satisfies Equation (2.41) can be used and 
an available check, as stated in [22], is that the eigenvector derivative is orthogonal to the 
eigenvector, hollowing, in trying to find an arbitrary non- trivial solution for|l/.}, the any 
kth row and column are partitioned from the eigenproblem matrix as follows: 
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(K-X,M\MxA fO' 

(K-^M)^^ r* r r 

Choosing k such that is non-zero and solving for and x^ in terms of Xf ,: 

\^ = -xA{K-A,M)i (2.50) 

Sinee the kth row (equal to the kth column due to symmetric matrices) is linearly 
dependent on the {n-1) remaining rows, the kth row can be eliminated and a system of {n- 
1) equations of rank (n-l) is guaranteed as follows: 

-(K-A,M\, (K-IMUM 

(K-X,M)„ (K-A,M\A\x-J ‘\{K-X,M\A 


Evoking that |l/.} is arbitrary; its kth element can be set to zero and using Equation (2.51) 
to solve for {if} results: 





A/f _ 

ADV ADV 

AK AAi 
ADV ADV 


M 

M 


Ml 

Ml 


1 


3 


(2.52) 


The unique solution for the eigenveetor derivative will then be: 



’'1 

< 0 > + c,4 

V ^ J 


(2.53) 
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where {V,} and {^ 3 } are calculated from Equation (2.52) and c. from Equation (2.48). 


An algorithm for calculating the eigenvector derivative (sensitivity), assuming eigenvalue 
sensitivity has already been calculated, is given below [ 22 ]; 

• Eoop over all eigenvalues 

• Calculate the left hand side matrix [D. ] = [ A - AM ] 


• Select k based on the largest element of } 

• Partition out kth row and column from [dJ 

• Decompose the matrix 1-^' ] 

• Eoop overall design variables 


M 


Calculate the RHS vector { F. | = - - 

^ ^ lADV ADV 

Eorward Backwards Substitution to calculate {if} 

Merge |i/.} to an n-vector by adding zero kth element 

Calculate c, =[m]{V,} 

Eigenvector derivative is given by } + c. 

• End of design variable loop 

End of eigenvalue loop 


F. CORRELATION OF STRAIN ENERGY TO SENSITIVITY MATRIX 


At this point, it would be useful to make a brief analysis of the strain energy of a 
dynamic system and to examine how it is correlated to the sensitivity matrix. A relevant 
analysis has been performed in [18], where the shape functions of an Euler- Bernoulli 
beam were used for the derivation of the strain energy of the system. In this work, after 
having set the EE model for the structure, which requires the calculation of the system’s 
mode shapes and of each element stiffness matrix, a simpler approach was made using 
the trivial equation; 
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U =-kx'^ 
2 


(2.54) 


which for each element of the model using modal data is translated into; 

=^[o"]^[a:"][o"] (2.55) 

where stands for the element stiffness matrix and stands for the partition of 

the mode shape matrix that eontains the eigenveetors related to the DOF assoeiated with 
the partieular element. It is of vital importanee to be noted here, that this mode shape 
matrix is not the original eigenveetor matrix of the baseline system, in the way that it has 
been modified in order to take into aeeount the physieal boundary eonditions of the 
system whieh may oeeur in some of the elements’ boundaries. For example, eonsider a 
eantilever beam with its left end elamped is examined and the eorresponding to the 

first element of the FE model. The entries of the matrix assoeiated with the right two 
DOF of the element will be eovered by the entries of the baseline system eigenveetor 
matrix assoeiated with the first two (uneonstrained as trivially defined) DOF of the 
eantilever beam, while the entries assoeiated with the left two DOF of the element must 
be fulfilled with zeros eorresponding to the eonstraint imposed to both the translational 
and the rotational DOF. 

That said, and having derived the stiffness sensitivity matrix for the same system, 
a very remarkable conclusion is reached. As shown in the following Figures 11 and 12, 
where all the modes of a 10 element eantilever beam are depleted, there is an exaet and 
direet eorrelation between the stiffness sensitivity matrix and the flexural properties of 
the strueture. An appropriate scaling with normalization to the maximum value of its 
mode has been made to both of them in order to make the eorrelation more legible. In 
these figures, the x-axis has the number of elements; the y-axis has the normalized value 
of both the sensitivity matrix and the strain energy, whieh are referred to eaeh mode of 
the system, respeetively. 
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Figure 12. 
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Correlation of Sensitivity Matrix to Strain Energy for Modes 11 Through 20 
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III. LINEAR SYSTEM SOLUTIONS IN SENSITIVITY BASED 

EEM UPDATING 


As demonstrated in the previous ehapters, the proeedure of sensitivity-based FEM 
updating involved the solution of a linear system, whieh was almost always 
underdetermined. This linear system is governed by the Equation (2.31): 

{Ary"}=[r]{ADE} 


In EEM updating, the veetor |ADy} eontains eorreetions to the model 
parameters, while in damage deteetion, the veetor {ADE} eontains flags to both the 
loeation and severity of struetural damage, relative to the undamaged EEM. In both 
oases. Equation (2.31) has to be solved for the design variables {ADE} . This process 
constitutes an inverse type of mathematical problem, which requires the inversion of the 
sensitivity matrix [T ]. Of course, this is feasible only if the matrix is square and full rank, 

which is rarely the case. Eor all cases where the sensitivity matrix is underdetermined, 
there are several algebraic methods of linear systems for solving the problem, which are 
covered in the literature [2]. 

A. USING BACKSLASH “\” AND “PINV” IN MATLAB 


In this work, Matlab is used to carry out this procedure. In dealing with these 
kinds of systems, there are two available commands/operators that can provide the 
required pseudo-inverse of the sensitivity matrix; these are “\” or “mldivide,” and “pinv.” 
Respectively, they use orthogonal-triangular (QR) decomposition, and singular value 
decomposition (SVD). Eventually, they both give a least-squares (square and full rank 
overdetermined cases) or a least-norm (underdetermined and rank deficient cases) 
solution to Equation (2.31), which is pursued by minimizing the norm of the residual, i.e.. 


min||[r]{ADE}-|A(y^} , except that they do it in different ways. Erom what has 
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already been discussed, it is obvious that for rank deficient problems, the solution 
obtained from each one is not exact. Some brief descriptions of how each of them is used 
[23] are given below. 

1. Using “\” or “mldivide 

The mldivide(A,B) and the equivalent A\B, perform matrix left division (back 
slash). A and B must be matrices that have the same number of rows, unless A is a scalar, 
in which case A\B performs element-wise division. If A is a square matrix, A\B is the 
same as inv(A)*B, where inv(A) stands for the inverse of A, except it is computed in a 
different way. If A is an n-by-n matrix and B is a column vector with n elements, or a 
matrix with several such columns, then X = A\B is the solution to the Equation AX = B 
computed by Gaussian elimination with partial pivoting. A warning message is displayed 
if A is badly scaled or nearly singular. If A is an m-by-n matrix and B is a column vector 
with m components, or a matrix with several such columns, then X = A\B is the solution 
in the least squares sense to the overdetermined system of Equations AX = B, only if the 
rank k is equal to n (number of columns). In other words, X minimizes norm(A*X - B), 
the length of the vector AX - B. The rank k of A is determined from the QR 
decomposition with column pivoting. If k < n or the system is underdetermined, this is 
not the same solution as x = pinv(A)*B, which returns a least-norm solution, but in this 
case X=A\B returns a solution that is called “basic” solution and has at most k non-zero 
elements per column. 

2. Using “pinv” 

The syntax B = pinv(A) returns the Moore-Penrose pseudoinverse of A. The 
Moore-Penrose pseudoinverse is a matrix B of the same dimensions as A' (transpose of 
A) satisfying the four following conditions; A*B*A = A, B*A*B = B, A*B and B*A are 
Hermitian. The computation is based on singular value decomposition of A. If A is 
square and not singular, then pinv(A) is an expensive way to compute inv(A). If A is not 
square, or is square and singular, then inv(A) does not exist. In these cases, pinv(A) has 
some, but not all of the properties of inv(A). If A has more rows than columns and is not 
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of full rank, then the overdetermined least squares problem of minimizing norm(A*x-b) 
does not have a unique solution. Two of the infinitely many solutions are x = pinv(A)*b 
and y = A\b. These two are distinguished by the faets that norm(x) is smaller than the 
norm of any other solution and that y has the fewest possible non-zero eomponents as 
described for the use of 

3. Example Demonstrating the Difference Between the Two Commands 

The following example is totally representative of the properties of the two 
Matlab operators for this kind of linear systems solving [23]. 

Assume the matrix generated by A = magic(8); A = A(:,l;6). This is an 8-by-6 
matrix that happens to have rank(A) = 3: 


64 

■ 2 

3 

61 

60 

6 


'260' 

9 

55 

54 

12 

13 

51 


260 

17 

47 

46 

20 

21 

43 


260 

40 

26 

27 

37 

36 

30 

, and the RHS is b = 260*ones(8,l), that is:h = < 

260 

32 

34 

35 

29 

28 

38 

260 

41 

23 

22 

44 

45 

19 


260 

49 

15 

14 

52 

53 

11 


260 

8 

58 

59 

5 

4 

62 J 


260 


The scale factor 260 is the 8-by-8 magic sum. With all eight columns, one 
solution to A*x = b would be a vector of all Is. With only six columns, the equations are 
still consistent, so a solution exists, but it is not all Is. Since the matrix is rank deficient, 
there are an infinite number of solutions. Two of them are x = pinv(A)*b, which is: 


x = t 


1.1538' 


'4' 

1.4615 


5 

1.3846 


0 

1.3846 

^, and y = A\b which produces the result: y = < 

0 

1.4615 


0 

1.1538 


-1 


with the following 


warning: Rank deficient, rank= 3, tol= 1.8829e-013. 
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Both of these are exaet solutions in the sense that nomi(A*x-b) and norm(A*y-b) 
are on the order of round-off error. The solution x is speeiai beeause norm(x) = 3.2817 is 
smaller than the norm of any other solution, ineluding norm(y) = 6.4807. On the other 
hand, the solution y is speeiai beeause it has only three non-zero eomponents. 

From the statement above, it is elear that in terms of damage deteetion and error 
loealization, the most preferable operator for solving for the linear system eorresponding 
to the strueture’s model is the “\”, as it returns the minimum number of non-zero 
eomponents. Using only one mode that is one sensitivity matrix row, it will return only 
one non-zero element, whieh for the right ehoiee of mode, will be the element with the 
true error (as will be shown later, the right mode is the one whose highest value happens 
to be at the element where the true error oeeurs); whereas “pinv” will return a veetor 
sealed by the mode used. Using more than one mode will provide the user with a more 
transparent result in terms of loealization until a square matrix is reaehed where there is 
no difference between the two commands, provided the matrix is not rank deficient. This 
may prove to be of vital importance during the pursuing of optimal configurations, 
especially because it makes clear that modes of systems/configurations with locally 
concentrated sensitivity values will reveal errors at the corresponding elements where 
these values are laid. 

B. SPATIALLY INCOMPLETE UNDERDETERMINED SYSTEMS 

Realizing from what is demonstrated above, the least square and/or the least norm 
sense in which Matlab approaches that inverse mathematical problem, it is essential that a 
thorough understanding of how this linear system solution is depicted to the dynamic 
system in terms of damage detection, is possessed. For this purpose, an FE model of a 
cantilever beam consisting of ten beam elements is used. First, the baseline model is set 
and the modal data of the system, i.e., natural frequencies and mode shapes, are 
calculated, along with the stiffness sensitivity matrix through a perturbation of one 
percent reduction to the rigidity of each beam element, recursively. Then, another FE 
model is built, identical to the baseline with the only difference being the presence of 
damage somewhere in the structure. This damage is represented by a 10 percent rigidity 
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reduction and this model is used to provide simulated modal data from a structure that 
has localized damage. We will refer to the damage as the "error," recognizing that these 
concepts apply equally to identification of model error, or structural damage. The main 
objective is, after having the data for both the baseline and the damaged model and using 
the stiffness sensitivity matrix, to examine how the model can predict the true element in 
error, which is the design variable that is subject to the 10 percent rigidity reduction. It is 
noted here that the procedure is repeated iteratively for all the beam elements, in order to 
assess possible errors anywhere in the model. 

1. Using Individual Modes 

The starting point of this analysis is the use of only one mode of the available 
data. For demonstration purposes, only the first five out of the total twenty modes will be 
taken into consideration. The rows of the stiffness sensitivity matrix [r] of the cantilever 

beam, which correspond to the modes of the system, are shown in the Figures 11 and 12 
of the previous chapter. The first five rows are reproduced in Figure 13. Again the values 
are normalized to the maximum value of each row: 
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Mode 1 



Figure 13. First 5 Rows of Stiffness Sensitivity Matrix 


The numerical representation of these first five modes is given below: 


0.0006 

0.0004 

0.0003 

0.0002 

0.0001 

0.0001 

0.0000 

0.0000 

0.0000 

0.0000 

0.0163 

0.0028 

0.0008 

0.0061 

0.0121 

0.0138 

0.0103 

0.0048 

0.0011 

0.0000 

0.0904 

0.0094 

0.0691 

0.0778 

0.0170 

0.0221 

0.0977 

0.1083 

0.0407 

0.0025 

0.2518 

0.1633 

0.2888 

0.0390 

0.2834 

0.2793 

0.0435 

0.3602 

0.3188 

0.0291 

0.5564 

0.7459 

0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 

0.1708 


Now, the first five natural frequencies of the baseline system and of the perturbed 
(damaged) system, and the first five entries of vector jAry^j formed by their difference 
are provided from the FE models. 

This is all the information needed for solving Equation (2.31) for the change in 
the design variables {ADR}. The following figures show the results for using a single 
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mode, i.e., individual modes 1 through 5, respectively, with the error being at each 
element, recursively. In all the figures, the y-axis represents the percentage error and the 
x-axis the number of the element: 















Actual Error 
Calculated Error 


10% Reduction at Elem # 1 


L 


4 5 6 7 8 

10% Reduction at Elem # 2 


L 


4 5 6 7 8 

10% Reduction at Elem # 3 


10 


10 



10 


Figure 14. Solution Using the First Mode Only 


35 













































































Actual Error 
Calculated Error 


10% Reduction at Elem # 1 
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23456789 10 

10% Reduction at Elem # 2 



1 23456789 10 

10% Reduction at Elem # 3 



Figure 15. Solution Using the Seeond Mode Only 
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Actual Error 
Calculated Error 


10% Reduction at Elem # 1 



10% Reduction at Elem # 2 




10% Reduction at Elem # 3 



Figure 16. Solution Using the Third Mode Only 
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Actual Error 
Calculated Error 



10% Reduction at Elem # 1 



1 23456789 10 

10% Reduction at Elem # 2 


5 S 7 9 iId 

10% Reduction at Elem # 3 



Figure 17. Solution Using the Fourth Mode Only 
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Actual Error 
Calculated Error 




10% Reduction at Elem # 1 



10% Reduction at Elem # 2 



10% Reduction at Elem # 3 



Figure 18. Solution Using the Fifth Mode Only 

Keeping in mind from Figure 13 that the maximum sensitivity values occur at 
elements 1, 1, 8, 8, 9 for modes 1, 2, 3, 4, 5, respectively, a very important conclusion is 
reached for the FE model solution using only individual modes, as this is clearly derived 
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from the study of Figures 14 through 18. The eonelusion here is that no matter where the 
aetual error is, the ealeulated result depends on and is exclusively led by the highest 
sensitivity value of the individual mode used for the solution. This covers the localization 
aspect of the error prediction. Moreover, it is observed that when the highest sensitivity 
value corresponds to the element being actually in error, then not only will the location of 
the damage be detected, but also, the true magnitude of the error will be revealed. This is 
the reason that there is an exact correlation of the calculated error to the actual error in 
cases like: first element in error using first or second mode, eighth element in error using 
third or fourth mode and ninth element in error using fifth mode. It can be predicted at 
this point that if, for example, only the nineteenth mode is used, which has the fifth 
element with the highest sensitivity value, only when the fifth element actually is in error 
will there be a positive result, which can be verified in Figure 19. 
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Calculated Error 


10% Reduction at Elem # 1 
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Figure 19. Solution Using the Nineteenth Mode Only 
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2 . 


Using More Than One Mode 


The next step of the current analysis will be the examination of the use of more 
than one mode, starting with the use of pairs of modes. Consider the first two modes and 
the corresponding stiffness sensitivity matrix rows from the two top plots of Figure 13. 
The numerical representation of them, denoted by [Tj ^ ] and using four significant digits 
is: 



0.0006 0.0004 0.0003 0.0002 0.0001 0.0001 0.0000 0.0000 0.0000 0.0000 

0.0163 0.0028 0.0008 0.0061 0.0121 0.0138 0.0103 0.0048 0.0011 0.0000 


while the natural frequencies of the baseline and the damaged system and their difference 
will again be calculated from the FE models. 

Now, solving for the design variables using the first two modes, we get the 
following figure, again looping over all the elements being in actual error recursively: 
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Actual Error 
Calculated Error 


10% Reduction at Elem # 1 


L 


4 5 6 7 8 

10% Reduction at Elem # 2 


10 



Figure 20. Solution Using Modes 1 and 2 

As shown in Figure 20, the solution always returns a vector with only two non¬ 
zero values (even with the last four elements being in error, there is a value at element 1, 

but it is relatively minor compared to that of the sixth element), since the rank of the 
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sensitivity matrix used, eonsisting of the first two rows of the full sensitivity matrix, is 2. 
Furthermore, these non-zero values always oecur at elements 1 and 6 of the beam. As 
explained above, it is no surprise that a positive result is obtained only when the actual 
error occurs at element 1 or element 6 or both, as will be shown later. In order to 
investigate the reason for this trend, the procedure for the solution will be examined in 
successive steps as follows: 

Consider again the stiffness sensitivity matrix of the first two modes shown 

above: 



0.0006 0.0004 0.0003 0.0002 0.0001 0.0001 0.0000 0.0000 0.0000 0.0000 

0.0163 0.0028 0.0008 0.0061 0.0121 0.0138 0.0103 0.0048 0.0011 0.0000 


This is a case of an underdetermined problem, which as discussed previously in 
this chapter, Matlab solves using QR decomposition with column pivoting. Column 
pivoting means that the columns’ norm is the factor taken into consideration for the 
pivoting and a permutation matrix is effectively applied to the columns of T. It is of vital 
importance to note here that higher modes have higher sensitivities and. Therefore, the 
higher modes tend to govern the norm of the sensitivity matrix column s . Keeping this in 
mind and for purposes of demonstration it can be considered, without any loss of 
generality, that a partial pivoting is adequate to explain how Matlab reaches this solution 
for this linear system. 


Moving on with this procedure, the largest entry that is considered the first pivot 
is entry 1 of the 2nd row. 



0.0006 0.0004 
0.0163 0.0028 


0.0003 0.0002 
0.0008 0.0061 


0.0001 0.0001 
0.0121 0.0138 


0.0000 0.0000 
0.0103 0.0048 


0.0000 0.0000 

0.0011 0.0000 


Performing a Gaussian elimination the matrix becomes: 



0 0.0003 
0.0163 0.0028 


0.0003 -0.0000 
0.0008 0.0061 


-0.0003 -0.0005 
0.0121 0.0138 


-0.0004 -0.0002 
0.0103 0.0048 


- 0.0000 - 0.0000 

0.0011 0.0000 
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The next pivot will be chosen from the first row, and it will be the entry with the 


highest absolute value, and it turns out to be entry 6; 



0 0.0003 
0.0163 0.0028 


0.0003 -0.0000 
0.0008 0.0061 


-0.0003 -0.0005 
0.0121 0.0138 


-0.0004 -0.0002 
0.0103 0.0048 


- 0.0000 - 0.0000 

0.0011 0.0000 


Performing the next elimination, the matrix becomes; 



0 0.0003 
0.0163 0.0130 


0.0003 -0.0000 
0.0091 0.0051 


-0.0003 -0.0005 
0.0019 0 


-0.0004 -0.0002 
-0.0006 -0.0004 


- 0.0000 - 0.0000 

- 0.0001 - 0.0000 


At this point, it can arbitrarily be said that the solution is governed qualitatively 
by the highest values of the final stage of the elimination of the rank-deficient matrix, as 
shown above. Comparing Figure 20 to this final stage, in terms of location, the solution 
for the design variables will always give non-zero values at the pivots, which also happen 
to have the highest values in each row. Again, remember that this is not absolutely the 
way Matlab manipulates its solution to this linear problem, but it serves to demonstrate 
the main ideas of this result and to illustrate how several model updating methods can be 
properly used. 


In order to validate what is stated above, modes 2 and 5 will be used together and 
an attempt to predict the result will be made. The distribution of these particular 
sensitivity rows over the beam is shown in Figure 13, while their arithmetic 
representation in a common matrix is given below; 



0.0163 0.0028 
0.5564 0.7459 


0.0008 0.0061 
0.2616 0.7508 


0.0121 0.0138 

0.5051 0.4994 


0.0103 0.0048 
0.7337 0.3113 


0.0011 0.0000 
1.1042 0.1708 


The first pivot will be entry 9 from the 2nd row; 



0.0163 

0.5564 


0.0028 

0.7459 


0.0008 

0.2616 


0.0061 

0.7508 


0.0121 

0.5051 


0.0138 

0.4994 


0.0103 

0.7337 


0.0048 

0.3113 


0.0011 

1.1042 


0.0000 

0.1708 
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From this point already, the first non-zero solution for the outcome should be expected at 
element 9. The elimination based on that pivot gives: 



0.0158 0.0021 0.0006 0.0053 0.0116 0.0133 0.0096 0.0045 0 -0.0001 

0.5564 0.7459 0.2616 0.7508 0.5051 0.4994 0.7337 0.3113 1.1042 0.1708 


The next pivot, that is the 1st row of the above matrix, will be entry 1: 



0.0158 0.0021 0.0006 0.0053 0.0116 0.0133 0.0096 0.0045 0 -0.0001 

0.5564 0.7459 0.2616 0.7508 0.5051 0.4994 0.7337 0.3113 1.1042 0.1708 


Even at this point, the second non-zero value of the outcome should be expected at the E* 
element, but for consistency reasons, the next elimination is presented: 



0.0158 0.0021 
0 0.6727 


0.0006 0.0053 
0.2422 0.5635 


0.0116 0.0133 
0.0965 0.0308 


0.0096 0.0045 
0.3955 0.1516 


0 - 0.0001 
1.1042 0.1751 


So, the prediction for this configuration is that no matter where the actual error is, 
Matlab will return a result with non-zero values at elements 1 and 9, which will be 
positive only when these elements are actually in error, either one or both. By the term 
positive result, the correct damage detection is meant in a qualitative aspect (error 
location), which might sometimes be consistent quantitatively (error magnitude). Yet the 
quantitative aspect cannot be proved, although in the cases analyzed in this subchapter, it 
is revealed. 

Solving now the fundamental Equation (2.31) for the design variables and looping 
again over all the elements that are in error recursively, we get: 
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Actual Error 
Calculated Error 


10% Reduction at Elem # 1 



10% Reduction at Elem # 2 



10% Reduction at Elem # 3 



10% Reduction at Elem # 4 



10% Reduction at Elem # 5 



Figure 21. Solution Using Modes 2 and 5 

It is easily shown in Figure 21 that what is stated above is validated exactly. 

In order to validate the presented procedure in a more generic manner, the first 
five modes will be used and again the result will be predicted without actually solving 
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Equation (2.31). A priori, the solution expected should have only 5 non-zero values, since 
the rank of the used sensitivity matrix, shown in Figure 13 is 5. In this case, all five plots 
of Figure 13 are used with the following arithmetic representation; 



'0.0006 

0.0004 

0.0003 

0.0002 

0.0001 

0.0001 

0.0000 

0.0000 

0.0000 

0.0000 


0.0163 

0.0028 

0.0008 

0.0061 

0.0121 

0.0138 

0.0103 

0.0048 

0.0011 

0.0000 

[7’.:5] = 

0.0904 

0.0094 

0.0691 

0.0778 

0.0170 

0.0221 

0.0977 

0.1083 

0.0407 

0.0025 


0.2518 

0.1633 

0.2888 

0.0390 

0.2834 

0.2793 

0.0435 

0.3602 

0.3188 

0.0291 


0.5564 

0.7459 

0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 

0.1708 


The first pivot 

in this case is entry 9 of the 5th row; 






'0.0006 

0.0004 

0.0003 

0.0002 

0.0001 

0.0001 

0.0000 

0.0000 

0.0000 

0.0000 


0.0163 

0.0028 

0.0008 

0.0061 

0.0121 

0.0138 

0.0103 

0.0048 

0.0011 

0.0000 


0.0904 

0.0094 

0.0691 

0.0778 

0.0170 

0.0221 

0.0977 

0.1083 

0.0407 

0.0025 


0.2518 

0.1633 

0.2888 

0.0390 

0.2834 

0.2793 

0.0435 

0.3602 

0.3188 

0.0291 


0.5564 

0.7459 

0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 

0.1708 


so one non-zero value should be expected at element 9. After the elimination based on the 
first pivot, the matrix becomes; 



'0.0006 

0.0004 

0.0003 

0.0002 

0.0001 

0.0001 

0.0000 

0.0000 

0 -0.0000 


0.0158 

0.0021 

0.0006 

0.0053 

0.0116 

0.0133 

0.0096 

0.0045 

0 -0.0001 

[7’.:3] = 

0.0699 

-0.0181 

0.0594 

0.0501 

-0.0017 

0.0037 

0.0707 

0.0968 

0 -0.0038 


0.0912 

-0.0521 

0.2133 

-0.1777 

0.1376 

0.1351 

-0.1683 

0.2704 

0 -0.0203 


0.5564 

0.7459 

0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 0.1708 

The second pivot will now be entry 8 of the 4th 

row; 





'0.0006 

0.0004 

0.0003 

0.0002 

0.0001 

0.0001 

0.0000 

0.0000 

0 -0.0000 


0.0158 

0.0021 

0.0006 

0.0053 

0.0116 

0.0133 

0.0096 

0.0045 

0 -0.0001 

[7’.:3] = 

0.0699 

-0.0181 

0.0594 

0.0501 

-0.0017 

0.0037 

0.0707 

0.0968 

0 -0.0038 


0.0912 

-0.0521 

0.2133 

-0.1777 

0.1376 

0.1351 

-0.1683 

0.2704 

0 -0.0203 


0.5564 

0.7459 

0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 0.1708 


so another non-zero value should be expected at element 8. After the elimination based 
on the second pivot, the matrix becomes; 
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0.0006 

0.0004 

0.0003 

0.0002 

0.0001 

0.0001 

0.0000 

0 

0 

0.0000 

0.0142 

0.0029 

-0.0030 

0.0083 

0.0093 

0.0110 

0.0124 

0 

0 

0.0002 

0.0373 

0.0005 

-0.0169 

0.1137 

-0.0509 

-0.0447 

0.1309 

0 

0 

0.0034 

0.0912 

-0.0521 

0.2133 

-0.1777 

0.1376 

0.1351 

-0.1683 

0.2704 

0 

-0.0203 

0.5564 

0.7459 

0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 

0.1708 


It must be said here that, in fact, there is no need to eliminate entries both above 
and below the pivots. The only elimination needed for the purpose of this demonstration 
is of the remaining rows not yet containing any pivot, which, as can easily be shown, 
yields to the form of an upper or lower triangular matrix after the appropriate column 
shifting. 

The third pivot will be entry 7 of the third row: 


fel= 


0.0006 

0.0004 

0.0003 

0.0002 

0.0001 

0.0001 

0.0000 

0 

0 

0.0000 

0.0142 

0.0029 

-0.0030 

0.0083 

0.0093 

0.0110 

0.0124 

0 

0 

0.0002 

0.0373 

0.0005 

-0.0169 

0.1137 

-0.0509 

-0.0447 

0.1309 

0 

0 

0.0034 

0.0912 

-0.0521 

0.2133 

-0.1777 

0.1376 

0.1351 

-0.1683 

0.2704 

0 

-0.0203 

0.5564 

0.7459 

0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 

0.1708 


making element 7 likely to have a non-zero value at the result. 

Again, after the elimination of the remaining rows without pivot, the matrix 
becomes: 


[ 7 ’.: 3 ] = 


0.0006 

0.0004 

0.0003 

0.0002 

0.0001 

0.0001 

0 

0 

0 -0.0000 

0.0107 

0.0029 

-0.0014 

-0.0025 

0.0141 

0.0153 

0 

0 

0 -0.0001 

0.0373 

0.0005 

-0.0169 

0.1137 

-0.0509 

-0.0447 

0.1309 

0 

0 0.0034 

0.0912 

-0.0521 

0.2133 

-0.1777 

0.1376 

0.1351 

-0.1683 

0.2704 

0 -0.0203 

0.5564 

0.7459 

0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 0.1708 


and now the next pivot is entry 6 of the second row, which makes element 6 one of those 
with a non-zero corresponding value: 
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'0.0006 

0.0004 0.0003 

0.0002 

0.0001 

0.0001 

0 


0 

0 -0.0000 


0.0107 

0.0029 -0.0014 

-0.0025 

0.0141 

0.0153 

0 


0 

0 -0.0001 

[7’.:3] = 

0.0373 

0.0005 -0.0169 

0.1137 

-0.0509 

-0.0447 

0.1309 


0 

0 0.0034 


0.0912 

-0.0521 0.2133 

-0.1777 

0.1376 

0.1351 

-0.1683 

0.2704 

0 -0.0203 


0.5564 

0.7459 0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 0.1708 

After the last elimination, that of the first row, the matrix becomes: 




'0.0005 

0.0004 0.0003 

0.0002 

0.0001 

0 

0 


0 

0 0.0000 


0.0107 

0.0029 -0.0014 

-0.0025 

0.0141 

0.0153 

0 


0 

0 -0.0001 


0.0373 

0.0005 -0.0169 

0.1137 

-0.0509 

-0.0447 

0.1309 


0 

0 0.0034 


0.0912 

-0.0521 0.2133 

-0.1777 

0.1376 

0.1351 

-0.1683 

0.2704 

0 -0.0203 


0.5564 

0.7459 0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 0.1708 

In this 

last matrix, the highest value of the first 

row occurs at entry 1: 




"0.0005 

0.0004 0.0003 

0.0002 

0.0001 

0 

0 


0 

0 0.0000 


0.0107 

0.0029 -0.0014 

-0.0025 

0.0141 

0.0153 

0 


0 

0 -0.0001 


0.0373 

0.0005 -0.0169 

0.1137 

-0.0509 

-0.0447 

0.1309 


0 

0 0.0034 


0.0912 

-0.0521 0.2133 

-0.1777 

0.1376 

0.1351 

-0.1683 

0.2704 

0 -0.0203 


0.5564 

0.7459 0.2616 

0.7508 

0.5051 

0.4994 

0.7337 

0.3113 

1.1042 0.1708 


Having reached the final stage of this procedure using the first five modes, the 
predicted solution is expected to have non-zero values at elements 1, 6, 7, 8, 9 and these 
values will constitute a positive result when the actual error occurs in each, or any 
combination of these elements. 

As far as the calculated solution from the solution of Equation (2.31) is 
concerned, the following figure is generated, which again validates what is stated above 
(again where the existence of non-zero values at elements 1, 6, 7, 8, 9 is not obvious due 
to the very low value relative to the other values): 
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Figure 22. Solution Using Modes 1 through 5 

It has to be recalled at this point that, although the quantitative aspect of damage 
detection (error magnitude) is present in the cases examined above, it is not proved in a 
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generic way for all situations. It will later be shown that there may be some discrepancies 
which distract the result from being positive as far as the error magnitude or the 
exclusivity of the location are concerned. These discrepancies will be dealt with in a later 
chapter. 


3. What if a Fully Determined Square Matrix is Available? 

From what is stated so far, a very important conclusion is reached. For a structure 
consisting of n elements, and considering one design variable for each element, if m 
modes are used for the sensitivity-based analysis, with m<n, then the maximum number 
of elements that can be positively identified to be in error is k, where k is the rank of the 
sensitivity matrix, which at most, is equal to m. In other words, if it is assumed that there 
is no linear dependency between the rows of the sensitivity matrix, then the maximum 
number of positive damage-detected elements is equal to the number of the modes used 
for the analysis. 

Therefore, it can be concluded that if damage detection throughout the structure is 
desired using one and only one sensitivity matrix, then this matrix must not only be 
square, but also full rank. This is the condition that the present work has been based on. It 
has to be said that even if a matrix is full rank, there may be some near-dependencies that 
are crucially negative for the desired result. This is why the enlargement of the data used 
must first be pursued and then the optimal configurations contained in this data must be 
found and used for globally positive structure identification. 

The following figures show that with a square matrix, comprised of either the first 
10 modes, or of the second 10 modes of the structure examined previously (cantilever 
beam), potentially good results are obtained (Figure 24) because all the elements may 
have a non-zero value in the solution vector. Apparently, this is only necessary, and yet 
not sufficient condition for positive damage detection (Figure 23), since different 
configurations give different results even though both are square and full rank. 
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Figure 23. Solution Using Modes 1 Through 10 
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Figure 24. Solution Using Modes 11 Through 20 
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IV. ENLARGEMENT OE TEST DATA- ENHANCEMENT OE 
SENSITIVITY MATRIX DETERMINABILITY/ CONDITIONING 


In the previous ehapter, the need for the eonstruetion of a fully determined 
sensitivity matrix for the deteetion of potential damage throughout the strueture was 
justified. In this ehapter, several methods will be pursued in this direetion until the point 
that optimal ABC eonfigurations are reaehed and prove more essential than other 
methods. 

A. USING EIGENVECTOR SENSITIVITY DATA 

In the theory ehapter, a method was presented for the ealeulation of the 
eigenveetor sensitivity matrix. This algorithm is applied to the same model that was 
analyzed in terms of eigenvalue sensitivity and the results are provided below. 

It should be mentioned at this point that even if data from eigenveetor are used, 
resulting in an improved FEM updating, two major drawbaeks as far as this analysis is 
eoneerned eannot be ignored. The first is the aeeuraey of the test measurements needed to 
aeeurately identify the mode shapes of a strueture sinee typieal measurement errors ean 
be eaused by noise, nonlinearities, ete. The seeond is that for large and eomplex 
struetures, the mode shape sensitivity analysis produees a large quantity of data, resulting 
in a more expensive problem in terms of ealeulation time and memory storage. 
Nevertheless, a possible eombination of the two aforementioned diserepaneies may eause 
the whole proeess to eonverge less rapidly or possibly eonverge into ineorreet results as 
deseribed in [24]. 

For this analysis, the same 10 element eantilever beam is eonsidered. Using the 

algorithm deseribed in Chapter II, a full eigenveetor sensitivity matrix is generated. The 

size of this matrix is 400x10. The number of eolumns is equal to the number of the design 

parameters, whieh again, are the rigidities eorresponding to eaeh one of the 10 beam 

elements. The number of rows eorresponds to 20 entries related to eaeh DOF of eaeh 

beam element, respeetively, whieh entries are generated for eaeh one of the 20 modes of 

the system. Obviously, the same FE model that was used for the derivation of the 
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eigenvalue sensitivity matrix is used again. The distribution of the eigenveetor sensitivity 
matrix for Mode #1, over all the elements, and for eaeh DOF is displayed in Figures 25 
and 26, while the first mode of the eigenvalue sensitivity matrix is available in the first 
plot of Figure 13. Eaeh plot is normalized to the maximum value of eaeh row for 
displaying reasons. It should be mentioned that in absolute values, the magnitude of the 
eigenveetor sensitivities is mueh less than that of the eigenvalue sensitivities, as shown in 
the arithmetie representation with four signifieant digits of the first mode of eaeh one 
below, and this may play a signifieant role in the proeedure of linear system solving: 

= 10 '^*[ 0.6034 0.4388 0.3017 0.1928 0.1117 0.0565 0.0234 0.0070 0.0011 0 . 0000 ] 


T =10“'’* 

eigvec,Model 


-0.2217 0.0842 
-0.0721 0.0274 
-0.5550 0.0334 
-0.0399 -0.0423 
-0.7109 -0.1599 
-0.0129 -0.0228 
-0.7201 -0.2476 
0.0090 -0.0071 
-0.6131 -0.2520 
0.0259 0.0050 
-0.4191 -0.1946 
0.0380 0.0136 
-0.1653 -0.0948 
0.0459 0.0192 
0.1247 0.0306 
0.0502 0.0222 
0.4322 0.1682 
0.0519 0.0234 
0.7449 0.3094 
0.0522 0.0236 


0.0588 0.0381 
0.0191 0.0124 
0.2235 0.1447 
0.0353 0.0228 
0.2373 0.3081 
-0.0285 0.0313 
0.0992 0.3243 
-0.0179 -0.0237 
0.0172 0.1980 
-0.0098 -0.0187 
-0.0235 0.0971 
-0.0041 -0.0152 
-0.0362 0.0132 
-0.0004 -0.0130 
-0.0320 -0.0605 
0.0016 -0.0118 
-0.0198 -0.1295 
0.0023 -0.0113 
-0.0052 -0.1971 
0.0025 -0.0112 


0.0224 0.0115 
0.0073 0.0037 
0.0850 0.0436 
0.0134 0.0068 
0.1806 0.0925 
0.0183 0.0093 
0.3024 0.1545 
0.0221 0.0112 
0.2960 0.2263 
-0.0220 0.0126 
0.1702 0.1985 
-0.0201 -0.0198 
0.0536 0.0816 
-0.0189 -0.0192 
-0.0577 -0.0330 
-0.0183 -0.0190 
-0.1665 -0.1465 
-0.0180 -0.0189 
-0.2746 -0.2597 
-0.0180 -0.0189 


0.0048 0.0015 
0.0016 0.0005 
0.0183 0.0055 
0.0029 0.0009 
0.0387 0.0117 
0.0039 0.0012 
0.0646 0.0195 
0.0047 0.0014 
0.0944 0.0285 
0.0052 0.0016 
0.1269 0.0382 
0.0056 0.0017 
0.0915 0.0483 
-0.0156 0.0017 
-0.0016 0.0198 
-0.0155 -0.0098 
-0.0943 -0.0391 
-0.0154 -0.0098 
-0.1870 -0.0979 
-0.0154 -0.0098 
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Figure 25. Mode 1 of Eigenvector Sensitivity Matrix for DOF 1 Through 10 
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Figure 26. Mode 1 of Eigenveetor Sensitivity Matrix for DOF 11 Through 20 
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The natural frequency and the eigenvector differences corresponding to the first 
mode can be derived from the baseline and the damaged model as before. In this case, the 
Model vector will have a dimension of 20x1. 

Two approaches will now be made to the identification of the effectiveness of the 
use of eigenvector sensitivity. The first one will use only mode shape data for the 
analysis, while the second will combine the same modes from both eigenvalue and 
eigenvector data. It is quite sufficient to use only the first mode for displaying the result 
for both cases as there are not many qualitative variations using other modes as well. 

Solving the fundamental Equation (2.31) for the design variables using only the 
eigenvector sensitivity data and looping over all the elements being in error recursively 
yields the following figure: 
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Figure 27. Solution Using Mode 1 of Full Eigenveetor Sensitivity Matrix Only 

It is pretty obvious that no positive result can be extracted from this approach, 
even though more than a square matrix, an over-determined matrix, is used for the 
analysis. 
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However, if eigenvalue sensitivity data is used along with mode shape data, a 
totally different outcome is produced. Using only the first mode, as seen in Figure 28, this 
configuration is capable of predicting damage anywhere in the structure: 
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Figure 28. Solution Using Mode 1 of Eigenvalue and Full Eigenvector Sensitivity Matrix 
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The same positive results yield if only the translational DOFs are considered to be 
measured and, thus, be considered in the mode shape data, as is shown in Figure 29. 
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Figure 29. Solution Using Mode 1 of Eigenvalue And All Translational DOF of 

Eigenvector Sensitivity Matrix 
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In spite of the positive results shown in Figures 28 and 29, it is realistically 
impossible for all DOF, or even all the translational DOF to be measured. Therefore, 
spatially incomplete mode shape data are examined in the following two ways. The first 
considers DOFs 3, 7, 11, 15, 19 of the cantilever beam to be measured and provides 
accurate mode shape data. The second one considers DOFs 3, 5, 7, 9, 11, 13, 15, 17, and 
19 to be measured. Then, the first mode of this data is used in each case for system 
identification through a combined eigenvalue and eigenvector sensitivity matrix, which 
are, respectively, a 6x10 underdetermined and a 10x10 square matrix. The solutions 
corresponding to these cases are demonstrated in Figures 30 and 31. The figures show 
that these are not configurations that give positive results no matter where the true error 
is. Actually, the first case positively identifies the error in six out of ten locations (which 
in terms of efficiency is the best it could do, since the sensitivity matrix is 6x10, as was 
described in Chapter III), and it will never predict the error in the remaining four 
locations, since the output vector has only six fixed non-zero components. On the other 
hand, although in the second case a square matrix is used, it is not able to predict the 
correct element being in error in all the locations. This is another example of why a 
square matrix is a necessary (first case. Figure 30) but not sufficient (second case. Figure 
31) condition for damage detection through a unique sensitivity matrix. 

Conclusively, it can be said that mode shape design sensitivity data may be used 
for the augmentation of the modal frequency data used for the sensitivity based FEM 
updating, while by themselves they do not constitute a reliable way to identify errors in 
the structure. Of course, considerable thought relative to the right choice of DOF to be 
measured should not be neglected, because different configurations have quite different 
effects on the output result as shown in the cases examined above. Beyond that, the 
drawbacks referred to in the beginning of this chapter must be taken into account in terms 
of the availability of measurement locations and of the accuracy associated with those 
measurements. 
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Figure 30. Solution Using Mode 1 of Eigenvalue and Translational DOFs [3, 7, 11, 15, 
19] of Eigenvector Sensitivity Matrix (Underdetermined) 
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Figure 31. Solution Using Mode 1 of Eigenvalue and Translational DOFs [3, 5, 7, 9, 11, 
13, 15, 17, 19] of Eigenvector Sensitivity Matrix (Square) 
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B. PARAMETER SELECTION 


As outlined in Chapter I under the background section, the first attempts to 
override discrepancies that are generated during the updating procedure were focused on 
what is called parameter selection. This method can be viewed in different variations, 
since several analyses have been made through the years and are provided in the 
literature. Some of them will briefly be described in the following paragraphs. 

1. Subset Selection 

The parameter subset selection method chooses a subset of the design variables 
that minimizes a penalty function based on the norm of the residual of Equation (2.31) as 
it is summarized in [5]. In that sense, it is assumed that only a subset of the parameters 
has corresponding independent columns of the sensitivity matrix and the optimum subset 
is chosen through an iterative procedure to produce a sub-optimal solution, commonly 
known as forward selection [3]. 

The first step is to find a single column among the columns of the sensitivity 
matrix that best represents the frequency differences vector Next, the combination 
of two columns that constitute the best sub-basis for the spanning of is determined, 
while the first parameter selected is retained. Then a Gram-Schmidt orthogonalization is 
performed and both the sensitivity columns and the Aty^are replaced by the 
corresponding vectors. The procedure is now repeated on this reduced problem for the 
next parameter that gives the smallest residual. An iterative process is then produced in 
order to combine the best parameters chosen with each one of those remaining to pursue 
the best data representation. It should be noted that the criterion of the smallest residual is 
equivalent to the smallest angle or distance of the vectors or the vector subspaces existing 
or created respectively. 

The subset selection problem may also be solved by using an exhaustive search of 
all possible subsets that minimize the penalty function and eventually choose the best 
subset that constitutes the global minimum of all the cases. However, taking into account 
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that the number of parameters and also the size of the subsets may be large for reasonable 
realistie problems, the eomputation required for an exhaustive seareh is prohibitive [5]. 

2. Clustering of Design Variables 

Another approaeh, similar to the subset seleetion method, is that dealing with 
elustering analysis. In this approaeh, the eolumns of the sensitivity matrix are eonsidered 
as veetors in multi-dimensional spaee. For example, if they are identieal, they appear as 
eo-linear veetors, meaning that the eorresponding parameters affeet the same way the 
modes of the strueture. The eigenvalue sensitivity matrix for a large FE model will 
eontain many quasi-dependent eolumns; these ean be identified via the elustering 
approaeh and may be eombined so that the eonditioning of the problem is improved. The 
most similar column s are added together, and in this way, the information available from 
the measurements is used to converge the model by the smallest adjustment to the 
updating parameters [9]. Again, in the general application of the clustering method, a 
distance criterion is needed to assess the similarity of different objects as is described in 
[10] and then as is also described, a method of linking the similar objects within the 
cluster is necessary. Clustering methods have been successfully demonstrated in [8], [9], 
[10], however, it has to be said that if during the updating, the sensitivities columns 
change, different clusters might be formed. Nevertheless, in the case of large and 
complex systems, clustering analysis is able to carry out the task of efficiently comparing 
many possible sensitivity columns. Together with sensitivity analysis and subset 
selection, it provides a method for including the effects of changes to less sensitive 
substructures by combining those columns of the sensitivity matrix that are close over the 
bandwidth of the modal testing. 

3. Grouped Elements Sensitivities 

A different manner in which selecting parameter is considered is presented in [7]. 
The concept of total sensitivity is introduced as the absolute sum of the sensitivities. The 
simple observation that this analysis is based on is the following. Suppose n design 
variables are considered and they are about to be merged into one that groups in a way 
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the effectiveness of the selected parameters into one sensitivity value. Then, if 
represents this value and represents the sensitivity value corresponding to each 
design variable, for / = 1; n, it is easy to see that: 

T,.r=T,+T,+- + T„=f^T, (4.1) 

/=1 

The total sensitivity, as defined above, will be: 

\tJ = \t,+t,+-+t„\<'£\t,\ (4.2) 

/-I 

The equality holds only when the signs of the individual sensitivity values of the 
design variables have the same sign. This leads to the conclusion that by grouping design 
variables into larger parameters, the number of parameters can be reduced, but the total 
sensitivity decreases, in general. From this basic observation, the parameter selection 
procedure, sets of parameters can be constructed with the objective that their effect on the 
residual of Equation (2.31) is maximized, meaning that the residual to be minimized is 
more sensitive to those selected parameters. 

This procedure is accomplished by a sequence of two different selection phases 
[7]. The first phase involves merging two neighboring elements into one parameter, as 
long as the signs of their sensitivity values are the same. Then, Equation (4.2) turns out to 
be equality and identical to Equation (4.1), so that the total sensitivity is just the sum of 
the merged parameters’ sensitivities. This phase ends when no neighboring parameters 
have the same sensitivity sign. 

If the number of the reduced parameters is not acceptable, then the second phase 
begins, but in this case, the sacrifice of total sensitivity is unavoidable. Obviously, none 
of the neighboring parameters have the same sign of sensitivities and by merging two 
neighboring design variables, the total sensitivity is changed: 

Zf.l+l^l+^hil^'.l+|r+ri| (4.3) 

k=l k=\ 

k^iyj k^i,j 

so that the decrement of the total sensitivity is expressed as: 
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( 4 . 4 ) 


\T\ + \T.\-\T +T. 

In other words, the sacrifice for reducing one parameter in terms of sensitivity, in 
this phase, is equal to Equation (4.4). This phase ends when some criteria that have been 
imposed, such as the number of design variables or the maximum allowed sacrifice, or 
both, are met. 

In this work, the concept of grouping neighboring elements is used to enlarge the 
available data. Since a cantilever beam is examined, and as seen in Figures 11 and 12, the 
sensitivity values are positive for all the design variables and for all the modes, no total 
sensitivity sacrifice is required. In other words, the first phase described above is 
sufficient for the analysis. So, the grouping of the 10 elements of the cantilever beam into 
5 pairs of neighboring elements results in the sensitivity matrix of dimension 20x5 whose 
distribution of all the rows is presented in Figures 32 and 33. Note, that in the following 
figures, the first bar corresponds to elements 1 and 2, the second bar corresponds to 3 and 
4 and so on until the fifth bar, which corresponds to elements 9 and 10. This means that, 
for example, the highest value of the first row of the grouped sensitivity can be taken to 
occur at either the first or the second element of a iO element cantilever beam. This will 
prove of crucial importance later on when the concept of grouped sensitivity will be used 
to enhance the use of ABC for this analysis. 
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Figure 32. Grouped Sensitivity Modes 1 Through 10 
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Figure 33. Grouped Sensitivity Modes 11 Through 20 
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c. 


THE USE OF ABC 


Before getting into the details of using ABC as an effieient way of enhaneing the 
sensitivity matrix determinability, it is useful to mention several things about the natural 
frequeneies of the ABC systems and their eonneetion to the baseline system natural 
frequeneies, with a short referenee to Cauehy’s Interlaee Theorem. 

1. Cauchy’s Interlace Theorem 

Interlace theorem, in general, deals with the eigenvalues of a principal submatrix 
related to those of the main matrix. According to the analysis made in [25], assume a 
matrix A of dimension nxn and the submatrix H of dimension mxm: 

with A’s eigenvalues: a^<a^<aj^---<a^ and H’s: 6^ <e,<e,- 

The theorem states that for y = 1,..., m it is: 

- ^nn-m and equivalently, < e_. < a_.. 

Alternatively, for = 1,2,..., n it is: 

^k-n.m and e_, < a_, < . 

In other words, d^j is an inner bound on relative to A’s spectrum. For example, 

‘0 £ o' 

and A= 8 0 1 Figure 34, 

0 1 0 

below, contains a schematic representation that demonstrates the theorem using the 
values of the above example. 


consider the following example, where H = 
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Figure 34. Schematic Representation of Cauchy’s Interlacing Theorem 

The actual application of this theorem to the ABC natural frequencies is that these 
frequencies will always interlace the natural frequencies of the system of higher order. 
This means that if a baseline system is considered, like the cantilever beam of the 
previous analyses, and one ABC is imposed on this system, then the eigenvalues of the 
ABC system will interlace the eigenvalues of the unconstrained cantilever beam. So, the 
first eigenvalue of the ABC system will be between the first and the second eigenvalue of 
the baseline system and so on. If a second ABC is imposed on the system, then the ABC 
system, with the two artificial constraints, will have natural frequencies that will interlace 
those of the ABC system with the one constraint. Therefore, in terms of bandwidth 
limitations, it is even better to use the ith mode of an ABC system than the ith+1 mode of 
the baseline system, if no other obvious reason for the use of each is present. Figure 35 
may be useful for the demonstration of this interlacing; 


♦ - ♦ - ♦ 


^^^ - Ei gen val u es of Basel i n e with 1 ABC 



♦ — ^ 


^^ Eigenvalues of Baseline 




Ascending Order 


Figure 35. Schematic Representation of Interlacing Between ABC and Baseline Natural 

Frequencies 




















2 . 


Using ABC 


Having said how ABC information can be obtained and where this information is 
relative to the speetrum of the baseline system, several ABC systems are used to faeilitate 
the error predietion in the eantilever beam that is examined in this work. All these 
eonfigurations are used in the eoneept of building an at least square sensitivity matrix that 
will be able to manipulate data from the FE model, and from the modal testing of the 
strueture, and provide a positive outeome via the Matlab linear system solution. 

The manner in whieh ABC are seleeted is based on the prineiple that sensitivity is 
inereased in areas near the artifieially imposed eonstraint, as shown in [18]. This 
observation ean help sinee the sensitivity distribution throughout the strueture is likely 
predietable, at least for the fundamental mode, of each ABC configuration. Therefore, the 
first step is to ereate 10 different ABC systems, whieh will consist of the baseline 
eantilever beam and a pin imposed on eaeh one of the 10 uneonstrained translational 
DOF, respeetively. All the data needed for these systems, like natural frequeneies, mode 
shapes and sensitivities, are obtained the same exaet way as for the baseline system. The 
first mode of the eorresponding sensitivities, along with that of the baseline system, is 
showed in Figure 36, where the red triangle depiets the DOF where the eonstraint is 
imposed and all the values are normalized to the maximum value of the first mode of 
eaeh system. It is easily noticeable how the eonstraints inerease the sensitivity near the 
DOF that they are imposed on, as long as these DOF are not near the free end of the 
eantilever. As one ean see,, the three last eonfigurations tend to simulate the seeond mode 
of the baseline system, whieh is not strange eonsidering the values of their natural 
frequeneies (elose to the seeond baseline eigenvalue) and their mode shapes (qualitatively 
elose to the second baseline mode shape). 
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Figure 36. First Mode of Sensitivity Matrix of Baseline (Plot 1) and of Eaoli ABC 

System (Plots 2-11) 


Next, the solution to the fundamental Equation (2.31) is pursued using two 


approaches. The first one uses all the configurations described in Eigure 36, thus 
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constructing an overdetermined 11x10 sensitivity matrix, while the second uses only the 
ABC configurations (Plots 2-11 of Figure 36), thus constructing a square matrix. These 
solutions are demonstrated in Figures 37 and 38, where the prediction of this model for 
all the elements being in error recursively is provided. 
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Figure 37. Solution Using Mode 1 of Baseline System and of Each ABC System 
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Figure 38. Solution Using Mode 1 of Each ABC System Only 

The two solutions provided in the above figures are quite similar since the 
sensitivity matrices used in each case do not have many differences in terms of sensitivity 
distribution along the structure, and since the frequencies difference vector used is almost 
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the same. What is also obvious, is that the use of these ABC eonfigurations is quite 
eneouraging sinee, at no eost, a relatively good solution is reaehed with the aetual error 
being predieted in all eases, both in loeation and in magnitude, the only drawbaek being 
the presenee of bias eonsiderably eomparable to the true error in some of the eases (when 
aetual element in error is 1,4, 5, 6, 7, 8). Towards the direetion of eliminating this bias, 
other eonfigurations will be pursued later in this thesis. 

A very important observation that is made from studying Figure 36 and imagining 
these plots as the rows of a sensitivity matrix is that, as was deseribed before, the last 
three rows seem almost identieal. This means that they add a degree of linear dependeney 
between them, whieh does not provide the linear system with any useful additional 
information, but also eonstitutes a diserepaney in the whole proeedure of the linear 
system solution as was stated in the introduetory ehapter. Therefore, only one of these 
three eonfigurations may essentially be used as a means of adding data to the system 
while the other two must be replaeed by other more useful eonfigurations. It ean be seen 
that even by exeluding the eonfiguration that eorresponds to a eonstraint imposed on the 
eighth translational DOF (plot 9 of Figure 36), whieh is the one closer to the second 
mode of the baseline system, and thus using a square sensitivity matrix, the solution of 
Figure 39 is achieved. This solution is somehow similar to those of Figures 37 and 38, 
but yet not best since there is still a degree of linear dependency between the last two 
rows. 
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Figure 39. Solution Using Mode 1 of Baseline and of Eaeh ABC System Exeluding the 
One That Constrains the Eigth Translational DOF 

Moreover, it can be observed that up to the point that the similar rows described 


above are excluded, the remaining sensitivity matrix presents a trend that reveals a very 
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large degree of linear independeney and well-conditioning and which is the milestone of 
the present work. This is the concept of a diagonally dominant or a tridiagonal dominant 
matrix that is described in the following subchapter and which will be sought through 
baseline, ABC and grouped element sensitivities configurations. This kind of matrix is a 
kind of banded matrix in terms of having its larger values on the main diagonal and on 
the first diagonal below and/or above the main diagonal. In sensitivity terms, it means 
that there exists a sensitivity concentration around all the different elements of the 
structure, uniquely identified in each row of the sensitivity matrix. This kind of matrix 
will be called diagonally dominant banded for the purpose of this thesis since, to the 
author’s knowledge; there is no official mathematical definition that can describe it. 

The first attempt at constructing such a sensitivity matrix using modes from the 
baseline system and from the ABC systems described so far is made using first’ fifth and 
twentieth modes of the baseline system (plots 1, 5 of Figure 11 and plot 10 of Figure 12) 
and the first mode of each ABC system that imposes a constraint to translational DOF 1 
through 7 respectively (plots 2 through 8 of Figure 36). That way, an effort to achieve 
high sensitivity values to all the elements of the cantilever beam in different sensitivity 
rows is made, along with keeping out rows that do not provide any useful information for 
the linear system solution, such as the three last plots of Figure 36. The sensitivity matrix 
then built is demonstrated in Figure 40 and the solution, with respect to the design 
variables, is provided in Figure 41. A second attempt uses the same configurations 
described in the first, replacing the fifth mode of the baseline with the second mode of the 
ABC system that constrains the eighth translational DOF. The sensitivity matrix 
constructed and the corresponding solutions are provided in Figures 42 and 43, 
respectively. 
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Figure 40. Sensitivity Matrix of Mode 1, 5, 20 of Baseline (Plots 1, 9, 10) and Mode 1 of 
Eaeh ABC System That Constrains Translational DOF 1 Through 7 (Plots 2-8) 
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Figure 41. 
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Solution Using Mode 1, 5, 20 of Baseline and Mode 1 of Each ABC System 
That Constrains Translational DOF 1 Through 7 
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Figure 42. Sensitivity Matrix of Mode 1, 20 of Baseline (Plots 1,10) And Mode 1 of 

Each ABC System That Constrains Translational DOF 1 through 7 (Plots 2-8) and 
Mode2 of ABC System That Constrains Translational DOF 8 (Plot 9) 
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Figure 43. Solution Using Mode 1, 20 of Baseline and Mode 1 of Eaeh ABC System 
That Constrains Translational DOF 1 Through 7 and Mode2 of ABC System That 

Constrains Translational DOF 8 


84 

































































As seen in the previous figures, the solutions provided by using these 
eonfigurations are quite aeeurate, both in loeation and in magnitude, as far as the error 
predietion is coneerned. Of eourse, these examples deal with a very simple strueture, for 
whieh several ABC eonfigurations with not only one eonstraint can be used, as long as 
natural frequencies of those systems are within the frequency spectrum available. 
Speaking of the available frequency spectrum, it should also be noted that modes higher 
than those implemented by the nature of each problem, which in general are a few lower 
modes, cannot be used as it is done with the twentieth mode of the baseline system in the 
above examples. Knowing this, and that the use of ABC systems may not be unlimited as 
mentioned, it is shown later in this work how other methods of enhancing the test data, 
like grouped element sensitivities, can be combined with the use of ABC systems, aiming 
at the construction of a diagonally dominant banded sensitivity matrix. 

3. Diagonally Dominant Matrix 

A short description of a diagonally dominant matrix is provided in [26]. 
According to that reference, a matrix is said to be diagonally dominant if in every row of 
the matrix, the magnitude of the diagonal entry in that row is larger than or equal to the 
sum of the magnitudes of all the other (non-diagonal) entries in that row. The concept can 
be extended to column entries. More precisely, the matrix A is diagonally dominant if: 

K' I - S K> I * ’ where denotes the entry in the ith row and jth column. 

It is very important to note here that if a linear system is examined where 
elementary row operations, such as switching rows, are allowed either through a 
permutation matrix or just by switching the corresponding entries of the right hand side 
vector, instead of examining the criterion of diagonal dominance for the diagonal entry, 
the maximum value of each row can be considered. Then, if the matrix under 
consideration has a maximum value at different entries for each row covering its whole 
length, with suitable row switching, it can take the form of a diagonally dominant matrix. 
Implementing this concept in the global sensitivity matrix used for model updating, if 
appropriate configurations for a structure can be found and combined together such that 
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in each row, different elements are subject to the highest sensitivity values while the rest 
of the values satisfy the criterion described above, this matrix can be said to be a 
diagonally dominant one. Furthermore, if the strict mathematical definition of diagonal 
dominance is not satisfied, but there is a significant concentration of high sensitivity 
values around the diagonal entries, then this matrix can be called a diagonally dominant 
banded matrix, as stated previously in this chapter. 

The main issue regarding this type of matrix is that its properties make it of 
particular use in several cases. It has been proved that convergence is guaranteed if a 
diagonally dominant matrix is used in many methods, like Jacobi or Gauss Seidel, for 
solving linear systems. Also, no partial pivoting is necessary for a strictly column 
diagonally dominant matrix when performing Gaussian elimination. By the Gershgorin 
circle theorem, a strictly (or irreducibly) diagonally dominant matrix is non-singular. This 
result is known as the Levy-Desplanques theorem [26]. These applications are of major 
importance in FE model updating since they guarantee that there is going to exist a 
solution to the linear system, which can be solved well-conditioned and with no concern 
about dependencies that may diverge the solution. The only concern is finding 
configurations that can sufficiently add data to the underdetermined problem, aiming at 
creating a diagonally dominant sensitivity matrix. 

D. MANIPULATING GROUPED ELEMENT SENSITIVITIES 

It has been mentioned above that grouping neighboring elements by simply 
adding their sensitivities values, with no sacrifice of the total sensitivity as long as they 
have the same sign, is an effective way of selecting parameters. Figures 32 and 33 show 
how the grouped element sensitivities are constructed for each mode of the cantilever 
beam examined, using pairs of two neighboring elements. The next thing to be 
considered is how these configurations can essentially be used to enlarge the available 
data and to assist in the effort of constructing a diagonally dominant banded matrix. 
Doing so can lead to a sensitivity matrix, which combines baseline, ABC and grouped 
element sensitivities information. 
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Consequently, it must be shown that grouped element sensitivities data can 
replace other data in the global sensitivity matrix. At first sight, this seems to be a 
difficult task, since from what is shown in Figures 32 and 33, if only one miscellaneous 
mode is examined, its length will be five, corresponding to the five pairs formed from the 
grouped elements. But the global sensitivity matrix where this grouped element 
sensitivity mode is to be used has length 10, corresponding to the number of individual 
elements of the structures which are not grouped. Note here that the desirable use of 
additional data is towards the enhancement of the determinability of the sensitivity 
matrix, and not towards the design variables increase, and that can be performed by using 
grouped element sensitivities as extra rows for the global sensitivity matrix. Therefore, 
there must be a way found for this five-component row to be permuted to a 10 component 
row and then added to the global sensitivity matrix, without any loss of its properties. 

For that purpose, the twentieth mode will be considered in for the concept used 
for that permutation to be explained. The highest value of this mode occurs in the last 
pair of the cantilever beam elements, which are elements 9 and 10. Assume that this 
value is denoted by ■ Then, the rate of change of the sum of the rigidities of 

elements 9 and 10, with respect to the change of the twentieth natural frequency will be 
equal torj.° 5 „„^io. In other words: 


rr<20 

gr,9andl0 


rj-<20 . rj-f20 

^el,9 ~^^eld0 


ADV, 


gr,9andl0 


Aco^ 


20 


ADV,, + ADV,,, 


Aco‘ 


20 


(4.5) 


From Equation (4.5) it can be seen that this particular sensitivity value can be 
used to characterize and identify the potential change in the design variable of element 9 
or 10 or both, due to the change of the twentieth natural frequency of the system. That 
leads to the conclusion that can be used as the sensitivity value corresponding to 

one of the elements 9 or 10, while the other will have zero sensitivity value, when both 
elements are considered in a sensitivity row. This can be equivalently expressed by the 
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fact that if, for example, one of the two rigidities of elements 9 and 10 is redueed by 10 
pereent, that is either ADI/,, or ADI/, jq is 0.1, then the total deerease to their group will 

be 10 pereent too. 

This conclusion allows a wide range of possible manipulations of the grouped 
element sensitivities data in order to be appropriately used as enhancement measures of 
the global sensitivity matrix determinability. Figure 44 demonstrates two possible ways 
of how the twentieth mode of the grouped element sensitivity can be used in the global 
sensitivity matrix (of length 10) for providing a sensitivity row with highest value at both 
element 9 and element 10, by simply adding zeros to the elements that are not required to 
have the sensitivity value desired. 

Permuted Grouped Element Sensitivity of Mode 20 with Highest Value at Element 9 



Permuted Grouped Element Sensitivity of Mode 20 with Highest Value at Element 10 



Figure 44. Mode 20 of Grouped Element Sensitivity Permuted to Attribute Flighest Value 
Either to Element 9 (Plot 1) or to Element 10 (Plot2) 
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E. 


COMBINING BASELINE, ABC AND GROUPED ELEMENT 
SENSITIVITIES FOR A DIAGONALLY DOMINANT BANDED MATRIX 


The ultimate purpose of this work is to eombine baseline, ABC and grouped 
element eonfigurations in a an essential manner into a global sensitivity matrix, whieh 
will be able to prediet ehange to the design variables of a system, given the differenee 
between the natural frequeneies of the system to measure (potentially damaged) and those 
of the original (undamaged) model. Having said that, the eonfigurations of Figure 40 will 
be used again, that is modes 1,5,20 of the baseline system and mode 1 of each ABC 
system that constrains the translational DOF 1 through 7, with two different alterations. 
In the first, mode 5 of baseline system will be replaced by the twentieth mode of grouped 
element sensitivities permuted to have the highest value at the ninth element (plot 1 of 
Figure 44), while in the second mode 20 of baseline system will be replaced by the 
twentieth mode of grouped element sensitivities permuted now to have the highest value 
at the tenth element (plot 2 of Figure 44). The concept of these replacements is to prove 
that modes of the baseline system that contribute to the diagonally dominant banded 
matrix with high sensitivity values at certain elements, can be mutually replaced by 
modes of grouped element sensitivities that have high values at the same elements, 
respectively. It is even better, under the principles instituted for the diagonally dominant 
banded matrix, if the replacing modes, apart from only high value at certain elements, 
have a significant sensitivity concentration around the element with the highest 
sensitivity value, but this is something that is achievable depending on the nature of each 
problem. 

Considering the first case mentioned above, the configurations used are presented 
in terms of sensitivity distribution per row in Figure 45, and the solution reached using 
these configurations is depicted in Figure 46. 
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Figure 45. Sensitivity Matrix of Mode 1, 20 of Baseline (Plots 1,10) and Mode 1 of Each 
ABC System That Constrains Translational DOF 1 Through 7 (Plots 2-8) and Mode 20 
of Grouped Element Sensitivities With Flighest Value at Element 9 (Plot 9) 


90 


























































































0 . 

- 0 . 

0 . 

- 0 . 

0 . 

- 0 . 

0 . 

- 0 . 

0 . 

- 0 . 

0 . 

- 0 . 

0 . 

- 0 . 

0 . 

- 0 . 

0 . 

- 0 . 

0 . 

- 0 . 


Actual Error 
Calculated Error 


10% Reduction at Elem # 1 




Figure 46. Solution Using Mode 1, 20 of Baseline and Mode 1 of Eaeh ABC System 
That Constrains Translational DOF 1 Through 7 and Mode 20 of Grouped Element 
Sensitivities With Highest Value At Element 9 
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It is obvious from Figure 46 that this combination of configurations works fine for 
predicting the error when it occurs in elements 1 through 8, but it is not positive for errors 
near the free end of the cantilever beam. It has to be mentioned here that even if 
sensitivity matrix is full rank, there is a discrepancy regarding the frequency difference 

vector . This discrepancy is related to the use of the same values in its last two 

entries, since the twentieth natural frequencies of the baseline and the experimental 
model are the same for both the original and the grouped element model. This is a very 
important observation because it leads to the conclusion that there is a limitation imposed 
when grouped element sensitivities configurations are used along with baseline ones. 
This limitation defines that the same modes of baseline and grouped element sensitivities 
cannot be used together in a global sensitivity matrix, because the subject of using 
repeated eigenvalues is raised and this produces discrepancies to the linear system 
solution. This limitation can be further applied to the use of ABC systems that have 
natural frequencies very close to those of the baseline system, as was discussed for the 
first eigenvalue of the ABC system that constrains the eighth translational DOF (plot 9 of 
Figure 36) in comparison with the second eigenvalue of the baseline system. 

Considering, now the second case where the twentieth mode of the baseline 
system is replaced by the twentieth mode of the grouped element sensitivities permuted 
to have highest sensitivity value at the tenth element. Figure 47 is generated by the 
sensitivity distribution along the rows of the sensitivity matrix created by these 
configurations. No repeated eigenvalues limitation is present and Figure 48 is produced 
to demonstrate the solution given from the linear system using this matrix. 
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Figure 47. Sensitivity Matrix of Mode 1, 5 of Baseline (Plots 1, 9) and Mode 1 of Eaeh 
ABC System That Constrains Translational DOF 1 through 7 (Plots 2-8) and Mode 20 of 
Grouped Element Sensitivities With Flighest Value at Element 10 (Plot 10) 
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Figure 48. Solution Using Mode 1, 5 of Baseline and Mode 1 of Each ABC System That 
Constrains Translational DOF 1 Through 7 and Mode 20 of Grouped Element 
Sensitivities With Highest Value at Element 10 
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Figure 48 provides a great result as far the loeation and the magnitude of the error 
is eoneemed, wherever in the strueture it may oeeur. It elearly demonstrates that baseline, 
ABC and grouped element sensitivities configurations can be essentially combined into 
one global sensitivity matrix for the structure, which has the properties of a diagonally 
dominant banded matrix as described previously in this work. With proper attention to all 
the limitations mentioned before, it is possible to reveal any potential damage in a 
structure or to update an FE model. 

In addition to what is concluded above, the same concept can be extended to cases 
where actual error occurs in more than one element of the structure. For this purpose. 
Figure 49 is generated and includes the solution using the configurations depicted in 
Figures 41, 47 and 45, respectively, with actual error occurring simultaneously in 
elements 1, 3 and 10. Again, it is proved that optimal configurations are those who adopt 
the concept of a diagonally dominant sensitivity matrix, but also considering issues 
generated by using similar modes of different systems. 
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Error using modes 1,5,20 of baseline and 1st of each 2-8 ABC syztem 



Error using modes 1,5 of baseline ,1 St of each 2-8 ABC syztem and 20th of grouped at 10th element 
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Solution For Actual Error in Elements 1,3, 10 Elsing: Plot 1-Configurations of 
41, Plot 2-Configurations of Eigure 47, Plot 3-Configurations of Eigure 45 
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V. CONCLUSIONS/RECOMMENDATIONS 


CONCLUSIONS 

1. The natural frequencies of the ABC systems can easily be identified from 
a square experimentally measured FRF matrix. These natural frequencies 
will interlace the natural frequencies of the baseline system in accordance 
with the Cauchy interlace theorem. The eigenvalue/eigenvector 
sensitivities of these ABC systems can be calculated using the same FE 
model as that used to generate the baseline system sensitivities. 

2. There is an exact and direct correlation between the stiffness sensitivity 
matrix and the flexural properties of the structure. This means that strain 
energy modes have the same distribution in the structure as their 
corresponding stiffness sensitivity modes. 

3. For overdetermined full rank problems, the use of Matlab operators “\” 
and “pinv” both return the same least square solution to the linear system 
created by the fundamental Equation (2.31), which minimizes the 
Euclidean norm of the residual, once with QR decomposition and once 
with the SVD method. 

4. Eor underdetermined or rank deficient problems, the solutions obtained 
from the use of those two operators are not the same. In this case, “pinv” 
gives a least norm solution, in the sense that the length of the solution 
vector is minimum, while “\” returns a basic solution with the significance 
of having a number of non-zero values equal to or less than the rank of the 
sensitivity matrix. 

5. Using only one mode for the identification of the structure, the “\” 
operator gives a solution that has one non-zero value corresponding to the 
element with the highest sensitivity value , while using the “pinv” operator 
returns a solution that is scaled by the sensitivity distribution of the mode 
used and eventually by the strain energy of the same mode. 

6. The use of “\” turns out to be more effective for system identification 
because it provides a more transparent result in terms of error localization 
and it makes clear that the use of modes with locally concentrated 
sensitivity values will reveal errors at the corresponding elements where 
errors exist. 

7. Solving underdetermined systems, that is using fewer modes than the 
number of elements potentially in error, will identify potential damage in, 
at most, as many elements as the number of modes used. This means that 
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an error can be uniquely identified anywhere in a structure only if a square 
full rank sensitivity matrix is used that is of a size equal to the number of 
elements potentially in error. 

8. Mode shape design sensitivity data may sometimes be used for the 
augmentation of the modal frequency data used for the sensitivity based 
FEM updating, as long as the availability and suitability of the 
measurement locations is appropriate, and the accuracy associated with 
those measurements is within the requirements tolerance. 

9. A useful means of compensating for an insufficient amount of test data is 
the grouped element sensitivities permuted in such a way as to match the 
sensitivity matrix length by simply adding zeros in elements that are not 
required to have a sensitivity value, and the total sensitivity of the group 
being assigned to the remaining member of that group. Grouped element 
sensitivities are not useful when the same modes of the baseline system 
are used, because issues of repeated eigenvalues and perhaps of mode 
quasi-dependencies are raised. 

10. Baseline, ABC and grouped element sensitivities configurations can be 
used together in one sensitivity matrix with the only limitation being the 
use of modes from different systems that might be close to each other, 
hence creating quasi-dependencies and close-to-repeated eigenvalues 
issues. 

11. If baseline, ABC and grouped element sensitivities are used to create a 
“global” square sensitivity matrix, which clearly concentrates high 
sensitivity values around every element in the structure in different modes, 
under the concept of, as described, a “diagonally dominant banded” 
matrix, then error can be identified anywhere in this structure both 
quantitatively, in terms of location, and qualitatively, in terms of 
magnitude. 

RECOMMENDATIONS 

1. Verily that the concept of the “diagonally dominant banded” matrix can be 
applied to the identification of errors associated with mass/density 
parameters. 

2. Create an algorithm/program that chooses the best combination from a 
variety of ABC and/or grouped element sensitivities sets, along with 
baseline sensitivities, into a sensitivity matrix that can identify a parameter 
error regardless of its location. This should be approached under the 
constraint of the availability of modes in a certain test bandwidth and the 
criterion of minimum computational time relative to the size and the 
complexity of the structure. If the conditions set above make this 
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procedure prohibitive, find eombinations of eonfigurations that create 
underdetermined sensitivity matrices, capable of damage prediction, each 
one in different areas/elements of the structure, as some regions of a 
strueture may be more neuralgie than others. 
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